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ABSTRACT: Perceptual studies suggest that the visual system uses the “rigidity” as¬ 
sumption to recover three-dimensional structure from motion. Tillman (1984) recently 
proposed a. computational scheme, the incremental rigidity scheme , which uses the rigid¬ 
ity assumption to recover the structure of rigid and nonrigid objects in motion. The 
scheme assumes the input to be discrete positions of elements in motion, under ortho¬ 
graphic projection. We present formulations of Tillman’s method that, use velocity infor¬ 
mation and perspective projection in the recovery of structure. Theoretical and computer 
analyses show that the velocity based formulations provide a rough estimate of structure 
quickly, but are not robust over an extended time period. The stable long term recovery 
of structure requires disparate views of moving objects. Our analysis raises interesting 
questions regarding the recovery of structure from motion in the human visual system. 
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1. Introduction 

An important source of three dimensional information is provided by the relative mo¬ 
tions of elements in the changing two dimensional image. The human visual system 
is capable of recovering structure' from motion, under both orthographic and perspec¬ 
tive projection, and in the absence of all other cues to 3 D structure (see, for example, 
Mile’s, 1931; Wallach ami O’Conne'il, 1953; llraunsteiu, 1976; hhanssem, 1978; Tillman, 
1979). In stuelying the- ceunputation e>f structure' fre>m motion:, erne' innne'diat.ely face's 
the problem] that the recemry e>f st ruct ure’ is unelercemst raineel; there' are' inlinite'ly many 
3 D structure's consist emt with a give>n pat.tern e>f motion in the- e-hanging 2 D image. 
Additional constraint is roepiiroel t,e» establish a uniepie' infeTpre’tation. 

Early pe’rce'pt.ual stuelie's sugge'ste'el that the' rigielity e»f objects may play a key re>le' in 
the' rerem'ry e>fstructure 1 fre>m motion (Wallae'h anel O’Cemne-ll, 1953; (libson ami Cibsem, 
1957; (Ire'em, 1961; .Johansson, 1975, 1977). Oompufaliemal stuelie's later e'stablishe'el that 
rigidity is a suliicie'ntly pewerful cemst-raint < e> ele'rive- a uniepie interpreted iem e>f struct ure, 
unele-r a variedy e>f viewing cemelitiems. For example 1 , Ullman and Frondiu (Oilman, 1979) 
shewed that uneler e>rthe)graphic project iem, three’ views of four nem coplanar points are 
suffkiont to guarantee' a uniepie 3 D inte’rpre’tatiem (uj) to an tuiavoielable re'fle'ctiem about 
the image plane). In the case of perspective projectiem, Longue't Higgins anel Prazelny 
(1981) proved that the instantaneous velocity lie-lei anel its first anel see-emel spatial deriva¬ 
tives at a point admit at most three different. 3 D inte'rpre'tations. Tsai and Huang (1981) 
shewed that, with the oxcoptiem of a few special eon figurations, two perspective views 
of seven pennts in motion are sufficient to guarantee a uniejue 3 D interpretation. Wax- 
man and Ullman (1984) also addressed the uniepieness of the recove'ry of structure under 
perspective projectiem, basing tiledr results on a kinematic analysis e>f continuems image 
flews. Adelitiemal the'ore'tical results have' been obtaineel for variems classes e>f restricted 
motiem, such as planar surfae-e's in motion (Hay, I960; Lemguet Higgins, 1984: Waxman 
anel Ullman, 1985; Ullman, 1985; Ne-gahelaripemr and Hern, 1985), pure translntory mo¬ 
tion (Clocksin, 1980), planar or fixe'el axis rotation (Hoffman and Fliuchbaugh, 1982; 
Webb and Aggarwal, 1981; Be>bick, 1983; Bennett and Hoffman, 1984a,b; Sugie anel Ina- 
gaki, 1984), and translatiem peri)enelicular to the rotation axis (Lemguet Higgins, 1983). 
A review of the the'ore'tical results regarding the> recovery e>f structure from motiem can 
be: found in Tillman (1983). 

From theere’tical stuelie’s of the structure from mediem problem, it can be eemclueleel 
that by exploiting a rigidity constraint, 3 D structure can he recovered from motion 
alemey, using image information that, is integrated ewer a small extent in space and in 
time. Tlie'se the'ore'tical stuelie’s have also give’ti rise to algorithms for deriving the rigid 
3 D structure of moving objects (for example, Ullman, 1979; Lemguet Higgins, 1981: Tsai 
and Huang, 1981). Experimentation with tlmse algorithms has reve-alcel two important 
limitations. First, although it is possible- in theory te> recover structure from medlem 
infennatiem that is iutegrate'd over a small extent in space' and time, such a strategy 
may ne>t be- robust in practice-. A small amount of error in the image- me-asurememts 
can lead to very differe-nt solutions (Ullman, 1983). Second, most provious algemthms 
elerive' a three dimemsionai structure only wliem a rigid iuterpre't.atiem is pejssible, and 
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otherwise do not yield any interpretation of structure or yield a solution that is incorrect, 
or unstable. 

The first observat ion above suggests that a robust algorit Inn for recovering structure 
should use motion information that is more extended in space or time. This conclusion 
is supported in recent computational studies by Nagahdaripour and Horn (1985) and 
tlllman (1984). Negahdaripour and Horn addressed the recovery of the motion of an 
observer relative t.o a stationary planar surface. It was shown that a robust recovery 
of both the observer motion and the orientation of the plane is possible when dense 
measurement's of the spatial and temporal derivatives of image brightness are integrated 
over a large 1 region of the changing image. Thus, consideration of motion information 
that is more extended in space can lead to a stable recovery of structure. Brass and 
Horn (1983) also proposed an algorithm for recovering the motion of an observer relative 
to a stationary scene, which integrates motion information over an extended region of 
the image. The study by Ullman (1984), which will be developed further in this paper, 
demonstrated that a robust recovery of st ructure is also possible when motion information 
is integrated over an extended period of time. The extension in time can be achieved, 
for example, by considering a large number of discrete frames or by observing continuous 
motion over a significant, temporal extent. 

With regard to the human visual system, the dependence of perceived structure on 
the spatial and temporal extent, of the viewed motion has not yet been studied systemati¬ 
cally, but the following informal observations have been made. Regarding spatial extent, 
two or three 1 points undergoing relative motion are sufficient to elicit a perception of 3 D 
structure (Borjesson and von Ilofsten, J973; Johansson, 1975), although theoretically the 
recovery of structure is less constrained for two points in motion, and perceptually the 
sensation of structure is weaker. An increase in the number of moving elements in view 
appears to have little affect, on the quality of perceived structure (for example, Petersik, 
1989). Regarding the temporal extent of viewed motion, Johansson (1975) showed that 
a brief observation of patterns of moving lights generated by human figures moving in 
the dark (commonly referred to as biological motion displays) can load to a perception of 
the 3 D motion and structure of the figures. Other perceptual studies indicate that the 
human visual system requires an extended time period to reach an accurate perception of 
3 D structure (Wallach and O’Connell, 1953; White and Mueser, 1900; Green, 1961). A 
brief observation of a moving pattern sometimes yields an impression of structure that is 
“flatter” than the true structure 1 of the moving object. (Tlllman, 1984). Thus, the human 
visual system is capable of deriving some sense of structure from motion information that 
is integrated over a small extent in space and time. An accurate perception of structure 
may, however, require 1 a more 1 e'xteneleel viewing period. 

It was noted earlier that most algorithms feu - recovering structure from me»tion are 
unable to interpret nonrigiel mediems. There 1 are 1 , howewe'r, some e'xce'ptiems to this. Bem- 
nett. anel 1 fedfmaii (1984b) stuelie'd the 1 minimum amemnt of motiou infortnatiem required 
to ele'rive a unique iut.erprotat.iou of the 1 structure 1 of a set e>f eliscrede 1 oleineuts undergoing 
nemrigid motiem, when it is assumeel that the elemiemfs are 1 rotating about a fixeel axis 
in space. Iloffman and Flinedibaugh (1983) prope>se>d an algorithm for interpreting the 
3 D motion anel structure in biological motiem displays. This algorithm decomposes the 




overall uomigid motion into pairs of points that, are rigidly linked and rotat ing in a plane, 
and triplets of points forming two hinged links that rotate' in the same plane'. Koenderink 
and Van Doom (198-1; Koenderink, 1981) examined the class of bending deformations, 
which satisfy the- physical constraint that distance's along the- surface' of the' object are 
pre'seTveel by the* t ranslbrmal iem. This e'lass of deformations excludes any strele-liing e>r 
e-emipressing e>f the e>l>je-e t surface. In its current, fornmlat iem, the- method proposed by 
Koenderink and Van I)e>e>rn for recovering the’ structure' e>f be'iieling surface's re'ejuires that 
the surface's lee- e-eunplete*, in contrast wit h other algorithms that are’ able- te» int-erpret the 
structure? of isolated points in motion. To cemcluele-, tlie’ algorithms elisentsse-el above for 
recovering the' structure of nonrigid objects in nie>tion all aehlress re'slriede-el e-Iasse’s e»f 
t.he'se motions, such as lixeel axis meitiem, planar meitiem, and be’iieling deformations. 

The mechanism for recovering structure fremi meitiem in the’ human visual system 
appe'ars neit. tc> he' base'el strie-tly on the' rigielity assumption. It is an e've'ryelay expe'rie'iice 
to pe're'e'ive' the structure anel meitiem of deforming objee'ts sue-h as a lleiwing rive>r, an 
expanding balloon, eir a elane-ing balle'rina. Perceptual stuelies reveal that the' human 
visual system can ele'rive' seune' se'nse* eif structure' for a breiael range’ of nonrigiel motions, 
including stretching, bemeling anel e've'n more’ complex type's eif eleformatiems (Johansson, 
J964, 1978; Janssem and Johansson, 1973; Teielel, 1982, 1984). Tt is also the case that 
displays eif rigiel objects in meitiem sennetimes give' rise' lei the' perception eif seunewhat. 
disteirting eibjects (Wallae'h, We-isz anel Aelams, 195G; White* anel Mueser, 19G0; Green, 
.1901; Braunstein, 19G2; Sperling et al., 1983; Hildreth, 1984a; Adclson, 1985). 

In this paper, we focus ein the recemt weak of Ullman (1984), which provides a 
morei flexible? metbeid for deriving the struct ure eif rigid anel nonrigiel objects in motion, 
and preivides a natural means for integrating motion informal ion over an extended time 
period. This method make's use of the rigielity assumption, but in a different, way from 
previous studies. The' algorithm, calleel the’ incremental rigidity schema , maintains an 
internal meidel of the structure of a moving object, which is continually updated as new 
positions of image elements are considered. The* initial mealed may be flat, if no either cues 
to 3 D structure are present, or it may be determined by other cues available, for example, 
freun binocular stereo, shading, texture', anel perspective (Marr, 1982; Ballard and Brown, 
1982; Horn, 1985). As oacli new view of the moving object appears, the algorithm 
computes a new set of 3-D coordinate’s for points on the object, which maximizes the 
rigielity in the tranfonnation from the current model to the new positions. In particular, 
the algorithm minimizes the change in the 3 D distance's between points in the model. 
The- formulation presented by Tillman assumes the input to the recovery process to consist 
of a sequence of discrete frames, each cont aining a set of discrete feature points. Through 
the process of repeatedly considering a new frame in the sequence and updating the 
current model of the structure of the moving features, tlie incremental rigidity scheme 
builds up and maintains a 3 D model, and can be applied to both rigid and nonrigid 
objects in motion. Further details of the incrementivl rigidity scheme arc presented in 
section 2 and in Appendix A. 

The incremental rigidity scheme has a number of advantages, from a computational 
perspective (Ullman, 1984): (1) because it integrates information over an extended time 
period, it provides a stable recovery of structure, particularly in the presence of error in 
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the- image- measurements. (2) il allows deviations from rigidity. while always maintaining 
some model of I lie 2 1) structure of Hie object, (.*5) it provides a natural means for 
interactions with other sources of ,‘5 1) information, and (1) empirical studies suggest, 
that the algorithm is able to recover the correct 3 1) structure. The behavior of this 
algorithm is also consistent in several ways with human perceptual behavior (l)llman, 
198-1). 

In this paper, wo develop Ullman’s work further, in several directions. First, in 
section 2, we present a continuous formulation of the incremental rigidity scheme, which 
uses velocity information at discrete points as input to the recovery process, in section 
3, we then examine in more detail, the behavior of the incremental rigidity scheme when 
presented with rigid objects undergoing rotation about a single axis in space, under 
orthographic projection. In particular, through computer simulations and a theoretical 
analysis, we examine the behavior of (he .discrete formulation as a function of the angular 
displacement between frames, and compare this behavior with that of the continuous 
formulation. Finally, in section 1, we present, both discrete and continuous formulations 
of the incremental rigidity scheme that use perspective projection. Through computer 
simulations, we begin to examine the behavior of the perspective formulations, when 
presented with rigid objects undergoing both pure' rotation about a single axis, and pure 
translation through space. 

The main conclusions of the paper are the following. The direct use of velocity 
information as input to the incremental rigidity scheme can provide a rough estimate of 
the structure- of a moving object over a slmrt, vie-wing period, but is not sufficiently perw- 
e-rfid to allow a ele-tailed ami robust rocem-ry of structure ewer an e-xte-nele-el time period. 
The computatiem of a stable long term sedation appe-ars to reepiiro the use erf views e>f 
a moving object that eliffe-r significantly. This implie-s the- nee-el for a recove-ry process 
with me-merry of the- past views, but this memory need not be extended indefinite-ly and 
continuously into the past. A small numbe-r of discrete views of a moving e>bje-ct are suf¬ 
ficient: for recovering 3 D structure. Tn the- case- e>f the incremental rigidity scheme, the 
use at every instant of a curre-nt, model of the- 3 D structure of the object, and a present 
view that is sufficiently different, from preceding views, can provide a robust recovery erf 
structure. In the case erf rot.at.iou erf a rigid object about a single axis in space, both the 
rate erf conve-rge-nce erf the- algerrithm ter the final sedation and the qualify erf the solution 
elecre-ase- as smaller angular elisplaceme-nts be-tween vie-weel frames are- cernsiele-re-d. In tlie 
limit erf the cerntinuerus formulation, the- solution is tier krngcr stable. The- behavior erf the 
perspective formulation of the incre-me-nfal rigidity scheme- is more complex than that 
erf the orthographic fornmlaliern. We- ferunel that, if the absolute persitiem of an erlrje-ct 
in space is knerwn throughout the motion of the erlrje-ct, them the* pe-rspectivo formula- 
tiern performs wedl, similar ter the- orfhergrajrhic ferrnndatiem. The- results erf computer 
simulations re-ve-aled a de-graelatiem in performance' with smalle-r angular and spatial dis¬ 
placements between framers, but this de-gradation was somewhat more se-ve-re than in the 
errthergraphic fonmdat.iem. Again, in the- limit erf the- continuous formidatiern, the serlut.iern 
is ner lemge-r stable-. Our analysis raises important questions regarding the quantitative 
ability with -which the- human visual system can recover structure from merfiern; those 
querstierns are discussed in section 5. 
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2. Discrete and Continuous Formulations of the Incremental Rigidity 
Scheme 

In this section, we first describe ITU mini’s (1981) original formulation of the incremental 
rigidity scheme', which assumes the visual input to consist of a. sequence of frames, each 
containing a number of discrete points that may correspond to identifiable' feature's in the’ 
changing image'. W<‘ then preseml a formulation that use's vedoeity informalie>n at eliscre-te 
points in a cemlinuemsly e-hanging image- as input, to the' re-eewe-ry pre>ce-ss. The- analysis 
in this section assume-s e>rthe>graphie- prefect iem of the 1 se-e-ne' emto the image- plane. 

The- me>t ivat.ions fe>r cemsiele-ring a e-emt inuems formulation are- thre-e-folel. First, on t lie- 
basis e>f the- re-sult.s e>f ee>mpute>r simulations, ITllman (198-1) ne>te-el that whe-n analyzing 
ohje-e-t.s unele-rgeiing rigiel rotatiem, the eemve-rge-nce of the- ine-re-me-nl al rigidity sehe-me to 
the- corre-e-t sedutiem was sleiwe-r whe-n smaller angular se-paratiems be-t.we-e'n frame's were- 
use'd. This suggests that the- se-he-me- may pe’rfeirm be-tte-r whe-n suce-e-ssive- vie-ws e>f the- 
ediject elilfe-r significantly. We- cemsiele-re-el the- limit e>f arbitrarily eleise- frame-s, both as 
a me-ans e>f stuelying this phe-nemienem, an el to ele-te-rmine- whe-the-r a robust re-cove-ry e>f 
structure- is still possible- unde-r tlie-se- e-emelitiems. A se-cemel motivatiem is that recent 
weirk e>u the- ceimput.atiem of an instantaneous 2 D ve-leie-ity fie-lel freun the- ehanging im¬ 
age- suggests that a uniejuo ve-leicity fie-lel can be- eibtaineel for ge-ne-ral classe-s e>f meition, 
exploiting a e-emstraint on the smoe>tlme-ss of the- ve-leicity fie-lel (Ileirn an el Sclnmck, 198.1; 
Hilelreth, 1984a,b; Nagel, 1984). Ultimate-ly, it may be- useful to integrate the results of 
such ve-le>city fie-lel computatiems with the; recovery e>f structure from motion. A third 
med.ivatiem is that III] man's formuhition of the incremental rigielity scheme leaels to the 
solutiem of a set of nemlinear equations. It is shown in Appendix A that, the continuous 
formulation presented he-re le-aels to the solution of a se-t of linear expiations. This makes 
a theoretical analysis e)f the sedutiem more accessible, ane’l could in principle result in a 
more; e-flicient compute-r impleme-ntatiem. 


2.1 Ullman’s Discrete Formulation 

The- incremental rigielity scheme maintains and upelat.es an internal moele'] M(t) of the 
vie-we-d object, which consists of a se't of 3-D coerelinati's: M(t) = (xi(t),yi(t), Zi(t)). 
In this se-etiem, we assume orthographie-, projection (the case of perspective- projection is 
addresseel in section 4) emto the X — Y image plane, so that (xi(t),yi(t,)) are the image; 
cemrdiuate's e>f the i th j>e>int, and Z{(t) is the current estimate- of the ele-pth at the i th 
pednt. (We assume a le-ft. handed coeirelinate system, with the positive- /, axis pointing 
away from the obse-rver.) In Ullman’s formulatiem, whe-n no other 3 D cue-s are present, 
the initial mealed M(t) at t — 0 is take-n to be flat; that is, Z{( 0) = 9 (or some other 
cemstant value) f<>r i — 1,... , n, whe-re n is the number of points in motiem. In principle, 
other initial cemfigurations coulel also be cemsiele-red. The the-ore-tical analysis e>f sertiem 3 
e-xamim-s the- lemg term stability of the ine-re-mental rigidity se-he-me, inele-pe-neh'nt of the 
initial moele-l of the structure- of the mewing points. 

Given a current moded M(t ) at time t, and the- image of the moving points in a 
new frame at a later time- t', the preddem is to e-emipute a new mealed M(t') suedi that 
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the transformation from A1 (l) In A l(t') is as rigid as possible. Since x L (l') and ?/,-(/.') are 
known, Ibis requires the computation of I lie unknown depth values -.,(/'). (It is assumed 
that tiie correspondence between points in the two successive frames is known.) The new 
depth values are computed as follows. Let lij(l) denote the distance between points i 
and j at time t. To make the transformation as rigid as possible, the values (/') for 
flit' new model are chosen so as to make /*_/(/) and lij(f') as similar as possible. For this 
purpose, Ullman defined a measure of the difference between hj{t) and /,-.,(//) as: 

= ill, (2.i) 

and formulated the recovery of structure as the computation of Zi(t') that minimize the 
following overall deviation from rigidity: 

D(U') = ^2d(l tJ (l,)J ij (t')). (2.2) 

by 

After the values Zi(t') have been determined using this minimization process, the now 
model M(t') (x,(//), yi(t'), [(■')) becomes the current model. A new frame is then 

registered and the process repeats itself. In this way, the scheme maintains rigidity by 
keeping the total distances between points in the model as constant as possible. The 
motivation for the cubic factor in the denominator of Eq. (2.1) is that; the nearest 
neighbors to a given point are more likely to belong to the same object than distant 
neighbors, so that a point is more likely to move rigidly with its nearest neighbors. The 
lf -(t) factor diminishes the influence of distant points in the recovery process. 

It should be noted that in the case of orthographic projection, only relative depth 
values, Zi{t) - Zj{t). can bo recovered, rather than absolute depth values, because under 
this form of projection, the image of a given object does not change with its absolute 
depth. In addition, 3 D structure is determined only up to a reflection about the image 
plane, since the orthographic projection of a rotating object, and its mirror image rotating 
in the opposite direction, coincide. 

2.2 The Continuous Formulation 

A continuous formulation, which uses velocity information at discrete feature points, 
can be developed as follows. Assume again that there always exists an internal model 
M[t) — (x;(f), yi{t), Zi(t,)). Assume also that the image velocities i,;(f) and are 

known. The problem is then formulated as the computation of the z components of 
velocity, z.;(t), that minimize the total continuous change in the distances between the 
points. Tin.' general form of the measure of overall deviation from rigidity is given by: 

D c (t)-^Y l d c (l ij (t)). (2.3) 

by 

In our analysis, we consider different possibilities for the measure, d c {lij(l)). The theo¬ 
retical development, of section 3.2, for example, considers the behavior of the incremental 
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rigidity scheme, as the frames become infinitesimally close, using the following discrete 
measure of t he change in the distance between pairs of points: 

d (/,,(/)./,,(/')) - (/ 2 (0 ( 2 - i ) 

The derivation of the continuous measure of rigidity is outlined in section 3.2.i. It results 
in the following expression for D c (t), as a function of the coordinates and velocities of 
the points (the arguments, t, have been omitted for simplicity): 


il(f) = ((*»• - *>)(*< ~ *>) +■ (vi - yjHm v-j) + (* *»•)(* - %)) 2 - ( 2 . 5 ) 

*>j 

This j>articu1ar measure of rigidity was used in the theoret.ical study for analytic sim¬ 
plicity. The comput<T simulations of the continuous formulation use the measure of 
change in the distance between points given again by the limit of (/,•;(<) - kj(?))*, lor 
infinitesimally close frames, which is: 

= (M*)) 2 . (2-G) 

In terms of the coordinates and velocities of the points, this yields the following overall 
measure of deviation from rigidity: 


D c (l) = 


x 3 -)(xj - Xj) T iju - yj)(yj - i/j) + {Zj - Zj)(zj 
{xi - xj ) 2 + {in - yj ) 2 + (Zi - Zj ) 2 


(Eqs. (2.5) and (2.7) use slightly different measures of rigidity, but serve the same role 
as measures of overall changes of rigidity for the continuous formulation. In the present 
paper we do not use different notations for these measures, as it will bo clear from the 
context which measure is used. More specifically, Eq. (2.5) is used in the theoretical 
analysis of section 3.2, and Eq. (2.7) is used in the computer simulations of section 3.1.) 
In other respects, the continuous formulation is similar to the discrete formulation. A 
model of the structure of the moving points is built up by continually taking into account 
new velocity information over an extended time period. Again, because orthographic 
projection is used, only relative velocities, 4(f) - Zj{t), can ho recovered. This can 
clearly be seen in Eqs. (2.5) and (2.7), in which the coordinate's of the points and their 
time derivatives all appear in difference's between pairs. Further details of the continuous 
formulation are presented in section 3 and Appendix A. 

The analyses presented in this paper mainly consider single rigid objects in motion, 
which are compact in the sense that the internal distances between pairs of points do not 
differ much from one another. In this case, the additional /.jh factor of Eq. (2.1) has little 
influence on the behavior of the' algorithm, so we omit f ed it in our theoretical analysis and 
computer simulations for the sake of simplicity. In general, however, a proper weighting 
(and not necessarily a cubic factor) of the influence of different distances among points 
is necessary for a belter performance of the algorithm. 
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3. Positions vs. Velocities as Input to the Recovery of Structure 

We stated earlier that on (lie basis of computer simulations, IJIlman (198-1) observed 
(bat when analyzing objects undergoing rigid rotation, the convergence of the incremen¬ 
tal rigidity scheme to the corn'd solution was slower when smaller angular displacements 
between frames were used. In this section, we analyze this phenomenon from a theoret¬ 
ical perspective^ focusing on the long term stability of the computed 3 D model. We 
first, examine the behavior of the continuous formulation, which uses a current model 
and measured image velocities at discrete points as input to the recovery process. We 
then turn to t he discrete formulation, which uses the discrete positions of discrete points 
as input, and examine the long term stability of its solut ion as a function of the angular 
displacement between frames. Our main conclusions are the following. First, for the 
particular class of motions considered, the discrete formulation always yields a 3 D so¬ 
lution that converges asymptotically to the correct solut ion, hut the rate of convergence 
varies with the angular displacement . The rate of convergence increases with increasing 
angular displacement up to a maximum, and then decreases with further increases in 
this displacement. The position of this maximum depends on such factors as the type of 
motion and geometric struct ure of the points. 

Although the orthographic projection is in general not physically valid, it is used here 
because' it allows a simpler formulation of the problem, and is therefore better suited to 
theoretical analysis and computer implementation. It, allowed us to gain a deeper insight 
into the nature of the phenomena studied. We nevertheless implemented the equivalent 
perspective formulation and confirmed that, the basic results remain valid under this 
formulation. The use of perspective projection enables the recovery of the structure 
of objects undergoing pun; translation, which was not possible under the orthographic 
projection. In the case of translation, the rate of convergence of the computed 3 D model 
to the true structure also increases with increasing spatial displacements between frames. 

Before presenting the results of the theoretical analysis, we illustrate the behavior 
of the discrete and continuous formulations through the results of computer simulations. 
We show that the continuous formulation yields an initial fast, convergence to a close 
approximation of the true structure of the moving points, but then oscillates over a large 
range. It consequently does not yield a stable long term recovery of structure. 

3.1 Observations from Computer Simulations 

In this section, we briefly illustrate' the behavior of the discrete' and continuous formula¬ 
tions of the ine-remeutal rigielity scheune, for the 1 special case of rigid relation of a small 
set, of discrete peeints about the vertical axis. (For eletails of the computer impleinonta- 
tion, see' Appenelix A.) In the case' of the eliscre'te feummlatiem, we' examine the rate' e>f 
convergence of the algorithm and the 1 ejuality e>f the solution that, it yie'lds, as a function 
e>f the angular displae-emenit be'twevn frame's. We- then e emipare it s pe'rformaucc with that 
e>f the' eemtinuous feermulatiem. In all *>f the examples pre'seuteel Iktc, the input, cemsisted 
of a set, of hve> points in spare. The- first- point is assume 1 *! to lie 1 at the eerigin *>f a <•<>- 
ordinate system that is displace’*! from the viewer along the 1 line erf sight,. The 1 position 
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of this first point is constant throughout the motion of the points. The coordinates of 
flic remaining four points were chosen randomly. Fig. I illustrates a. typical set. of five 
points, showing (heir project.ions onto the A' Y plane (Fig. In) and the A" Z plane 
(Fig. lh). In the simulations, the project(*<1 positions and velocities of the points wore 
computed analytically, rather than measured from real image sequences. 




Figure 1. A set of five points with random coordinates is projected onto (a) the X — Y plane, 
and (b) the X — Z plane. 


Fig. 2 illustrates the behavior of the discrete formulation of t he incremental rigidity 
scheme, as a function of the angular displacement between frames. Each figure shows 
a birds’ eye view of the set of rotating points (that is, their projection onto the X — Z 
plane), with filled circles representing the true positions of the points and open circles 
representing the structure computed by the algorithm. Fig. 2a shows the true positions 
of the points and the initial model at time t = 0. The initial model is assumed to be 
flat. Figs. 2b and 2c show the true and computed configurations of points after .120° of 
rotation. This final position was reached by taking three steps of 40° (Fig. 2b), and ,12 
slops of 10 (Fig. 2c). It can be seen that the use of a smaller number of more disparate 
views yields a 3 D model that is closer to the true solution. As noted in Appendix A, a 
steepest descent minimization algorithm was used for most of our computer simulations. 
We also analyzed a small number of examples using an exhaustive search algorithm to 
find the minimum solution, and again found the accuracy of the solid ion to vary with 
the angular displacement, between frames. We therefore believe this to he a fundamental 
behavior of the incremental rigidity scheme that is not simply a consequence of the 
particular algorit hm used to implement the scheme. This observation is supported further 
by the theoretical analysis of the next section. 

In Fig. 3, we show a series of graphs that illustrate in a different way, the behavior 
of the discrete formulation of the incremental rigidity scheme as a function of angular 
displacement. In this ease, the set of points shown in Figs. L and 2 were rotated by 







four full revolutions. Each of the graphs show the error between the true and computed 
structures, as a function of time. In particular, the following quantity is plotted: 

1 
2 

(3.1) 

where dij is the correct, 3 D distance between points i and j in the object, and l,[j is 


y>i - in ) 2 






Number of Revolutions 







Figure 3. Graphs of the error.in the internal distances between points in the computed 3 D 
model, as a function of time. The points are rotated 4 full revolutions, in steps of (a) 40", (b) 
20 '. (c) 10°, (d) 5° and (e) 1°. (f) The graphs in 3a 3e are shown superimposed. 


1 
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llit' corresponding distance between points i and j in tlit' computed model. The graphs 
shown in Figs. 3a through 3e correspond to rotations with angular displacements of 40", 
20', 10", 5 and 1", respectively. In Fig. 31, the live graphs are shown superimposed. 
Again, it. can he se(m that the rate of convergence and quality of t he solut ion improves 
with larger angular displacements. For a particular total angular extent, the error in the 
computed model decreases with increasing angular displacement bet ween frames. In the 
case of 40" displacements, the algorithm converges asymptotically and monotonically 
to the final solution. For smaller displacements, the convergence is no longer strictly 
monotonic, but is still essentially asymptotic toward the final solution. Fig. 4 shows 
the same set of graphs superimposed, but with !lit' error plotted on a log scale. The 
convergence of the solution is now essentially linear, with varying slopes, suggesting that 
the actual convergence is exponential. 



Figure 4. Graphs of the error in the internal distances between points in the computed 3 D 
model, as a function of time, plotted on a log scale. The points are rotated 4 full revolutions, 
in steps of 40°, 20°, 10°, 5° and 1°. The graphs are shown superimposed, with the angular 
displacements indicated above each graph. 


Fig. 5 illustrates the behavior of the continuous formulation of the incremental 
rigidity scheme. The same set of points used previously was again rotated about the 
vertical axis and the 3-D model was computed at infinitesimally closely spaced times, 
using the instantaneous velocities projected onto the image (sec' Appendix A for details). 
In Fig. 5a, we compare the true positions of the points (filled circles) with the best 
solution (open circles) obtained over 10 full revolutions of the points. Although the model 
is quite close to the true structure at this position of the points, the solution oscillates 
significantly over an extended time period. A graph of the error in the computed model 
over the' 10 revolutions is shown in Fig. 5b. There is an initial fast convergence toward 
the true structure of the points, but the algorithm then oscillates with high amplitude 
away from the true solution. When the error in the model is high, depth reversals often 
occur. Fig. 5c shows an example of the true and computed structures at a time of 
complete depth reversal. Such reversals were also observed with the discrete formulation 
of the scheme, although rarely. 

Fig. G illustrates the typical behavior of the discrete formulation of the incremental 
rigidity scheme, averaged over 10 configurations of five points. For each of the eonfigura- 
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Figure 5. (a) The best solution (open circles) obtained over 10 revolutions of the points is 

compared with the true positions (filled circles), for the case of the continuous formulation, (b) 
The error in the internal 3 D distances between points in the model, as a function of time, (c) 
A complete depth reversal between the true and computed structures. 


tions, Hie first point was placed at the origin of (he coordinate system and the coordinates 
of the remaining four points were chosen randomly. The set of points was then rotated 
by three discrete angular steps, with the size of the angular displacements ranging from 
1° to 90”. The initial 5 1) model for the points was Hat, and a new model was computed 
for each of the three discrete positions of the points. After the third step, we computed 
tli<> following measure of the absolute error in the internal distances between points: 
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1 \ ' I dij hj] 

n 2— dij 
* -.) 

wlicrc dij is 1.1 io (mo 3 1.) distance between points i and j, l rj is (ho distance between 
points i and j in the computed modol and n is tho total ninnhor of pairs of points. 
This moasnro was olioson boranso it expresses an average of the orror in tho modol 
relative to (ho t.ruo structure. (Nolo that a lower valno for this measure corresponds to 
loss (Ti’or in tho oompntod modol.) W<' also considered other measures of c>rror in tho 
distances between points and in tho actual depth values, and found the general behavior 
of tho algorithm to be tho same under different measures. Tho above' error moasnro was 
averaged over the JO configurations of points, and is plotted in Fig. 0, as a function of 
angular displacement. it can iirst be soon that in general, tho orror after three discrete 
slops of the algorithm varies with the size of the angular displacement between frames. 
There is a steady improvement in performance as tho displacement increases, to about 
50", followed by a degrcdal ion for increment s of GO", and a steady decrease in performance 
from 70" to 90". Tho degradation in performance for an angular displacement, of GO" was 
common; 8 of the 10 configurations of points exhibited t his behavior. This degradation 
may be a consequence of the symmetry between the initial and final views, which are 
rotated 180" from one another. In general, the convergence of the algorithm for three 
discrete steps degrades significantly for smaller angular displacements. This result is not 
surprizing, in that there is very little change in the discrete views for such small angles 
and a reduced total angular extent. The deterioration for large angles probably occurs 
because at 90", the number of views available to the structure from motion process is 
reduced to two. 



Angular Displacement 

Figure 0. Absolute error in the internal distances between points in the 3 1) model computed 
by the incremental rigidity scheme. The graph shows the error after three discrete views of the 
rotating points are considered, averaged over it) random configurations of five points. 
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Fig. 7 illuslrail's I,lie error in (lie computed model, as a function of angular dis¬ 
placement, for (lie case in which the overall rotation of the points was kept essentially 
constant . The same set, of 10 random configurations of points was rotated by discrete 
angular stops, with the size* of the steps varying from 10" to 00 . In this case, each 
configuration was rotated by a total of (approximately) 180' and 3G0 , and the same 
measure shown in Eq. (3.2) was then computed. This measure' was again averaged over 
the 10 configurations of points, and is plotted as a function of angular displacement, in 
Fig. 7. Fig. 7a shows the data for the case of 180 of rotation. Note that for the angles 
40", 50”, 70" and 80", there are two points plotted, corresponding to multiples of the 
angular displacement that are'just le\ss than and greviter than 180 . The' graph that, is 
supe'rimpe>se*el on the’ points passe's between the' pairs erf points for tlie-sc angle's. Fig. 7b 
shews the> same 1 elat.a for the’ case' of .‘>00 of rotation (fe>r the' ease 1 of angular displacements 
erf 50 1 anel 70", the' points were rot at e>el by a total e>f 350"'). Alsei shewn in Figs. 7a anel 
7b are- the average- e-rrors in the sedufiem that we-re eleriveel by the continuous formulation 
erf the> algeirithm aftor 180 anel 300 e>f rotatiern of the- ft) configurations of points. These 
two data perints are- inelicate'el by the stars along the' orelinate' of the graph. The- main 
obse-rvatiem to be made is that when the- points are- rotated by a earn slant total amount, 
there' is still a stremg depenelence erf the' rate’ erf cernvorge’noo erf the solution on the size 
erf the angular displacennent be-tween frame's. Tlie-re is again an imprervoinont, in conver- 
ge-nce rate as this angle incre'ases, up ter about, 50°, folleweel by a slight worsening ferr 
an angle' of 60°, then improveine’nt, again ferr 70", follewe'd by a steady werrsening to 90°. 
The elegraelat.ion in performance ferr an angular displacennent erf G0 r was again common, 
oce'uring ferr 7 erf the 10 conliguratiems of points. The ele-terierration of the convergence 
rate with decre>asiug angular elisplacemonts is not erlrvions, because in spite of the fact 
that the> change's between consecutive frames are smaller, tlie-re 1 are many more frames 
altogether. In aelelitiem, while the continuous ferrmulatiern provides a gererel estimate of the 
true solution after 180" (see Fig. 7a), the solution then degrades significantly, providing 
a relatively poerr serlutiern after 300° erf rertatieru (se-e Fig. 7b). 

Ter cerncluele', the 1 re'sults erf computer experiments with the eliscrete formulation erf 
the incremental rigidity se'hetne shew a e-le-ar depemle’uce of the' behavior erf the’ computed 
solution on the size- erf the angular displacements bet.weeu frames. In the limiting case 
of the continuous ferrmulatiern. there is an initial fast convergence of the solution ter war el 
the> true serlutiern, follerwerl by a substantial oscillation of the solution. The next se'ction 
presents a theore'tieal analysis that attempts ter explain this phenomenon. 

3.2 Analytic Study of Convergence Properties 

In this sc'e'tiern we- emtlinc the main eemelusious erf our theoretical analysis. The purpose erf 
this analysis is not a general study erf the behavior erf the increrne'ntal rigielity scheme, as 
sueh an analysis would be- too euiinherserme. Rather it cerne-ent rates ern a ferrmal analysis 
of the> cernvergence proport ie-s erf the- algorithm, ferr a family erf example's that we> he-lievve 
are- relevant, to the' ge'neral re'cervery erf struct lire- frerm motion. We emphasize the* roneepts 
raise-el by this analysis, anel cernune-nt ern ilu'ir ge'nerality. This se'ction cernsiders a subset 
erf those 1 example's analyzeel in the cermpute>r simulations of se’ction 3.1. 
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Angular Displacement 

Figure 7. Absolute error in the internal distances between points in the 3 D model computed 
by the incremental rigidity scheme. The graph shows the error after a total of (a) 180° of 
rotation, and (b) 300° of rotation. The errors are again averaged over 10 random configurations 
of five points. The stars along the ordinate indicate the data points for the average errors of 
the continuous formulation after 180° and 300°. 

To simplify the theoretical analysis, we consider in this section a measure of overall 
deviation from rigidity that is somewhat different from the one used by Tillman, and it 
is the following: 

Ai(M') = XXW - l W) 2 - (3.3) 

hi 

This new measure uses the difference between the squares of the distances, /?•(/), rather 
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Mum Uk- dilferemce lx-twe-e-n I.lie /,y(f)'s llie-mselve-s, Mins avoiding Mio use e>f sepiaix- roots. 
This measure also does not contain the cubic factor, /; ! ; (/). in the denominator, which was 
included in the measure proposed by Mllman as a means of reducing the contributions 
of distant points relative to nearby ones. In this theoretical analysis, we mainly consider 
configurat ions of points that are compact, in t he sense that the internal distances bet ween 
pairs of points do not. differ much from one another. In this case, the cubic factor seems 
not to be important. 

The main problem that this theoretical analysis addresses is the long term stability 
of the solut ions obtained by the incremental rigidity scheme, as a function of the angular 
separation between successive 1 frame's, for the- case 1 e>f rotation e>f rigiel objects about a 
single- axis in space. In the- pre-viems se-ctiem, we- obse-rve-el the- variatiem of e-emve-rge-nce 
rate with angular elisplaeenu'nt, in the* re-sults e>f compute-r simulatiems. In the- pre-se-nt 
aualysis, we- examine lliis e-on verge-nee- phe-uemu-nem through a stuely e»f the 1 stability of 
the incre-me-utal rigielity sche-me-. in particular, we- analyze 1 the- behavior e>f the 1 internal 
mexle-I unele-r small pe-rturbations, wlum it is ne-ai - the true 1 soluMem. We- first e-xamine the 
limiting eiase, whe-re the frames are- arbitrarily close- te> one- anedhe-r, that is when V t. 
The 1 analysis e>f this e-ase is pre-sente-el in se-ctiems 3.2.1 an el 3.2.2. We 1 then explore the 
cemve-rge-nce prope-rties e>f the 1 eliscre-t.e formulation in sect.iems 3.2.3 anel 3.2.4. Finally, in 
se-ctiem 3.2.5 we- summarize- the- Mie-ece-tie-al preipe-rtios eif the- two scheme-s, as a fune-tiem 
e>f angular elisplace-me-nt. 

3.2.1 The Continuous Formulation 

In the- case of the continuous fecundation, can be- ele-teirmine-d to a good first eirelc-r 

approximatiem by /,;/(£) and its derivative-, lij(t), that is: 

+ x (f -t). (3.4) 

As describe-el e-arlie-r, when such an approximatiem is used, rather than computing the- 
ele-pth values Zi[t') that maximize rigielity, we earn compute the- z coni{)ements of velocity 
Zi{ t) that elo the- same task. It earn be shown that in the- limit, as t' —► t, using Eq. (3.3), 
anel the eh-finitiem of l tJ , the qufintity to be- miuimize-el is: 

Dc{t) * ~ *j)(**‘ “ + (?h ~ V]Min ~ Vi) + i z i ~ z j)i z i ~ z j )) 2 - (3-5) 

Thus, giveai the mode-1 at time t, and the x anel y compone-nts of velocity in the image, 
we compute the- temporal ele-rivative-s of the ele-pth values that, minimize D c (t) give-n by 
Eej. (3.5). This is Mie- cemtinuous formulatiem that we use- in the theoretical analysis of 
the convergence preipe-rtios of the: incremental rigielity scheme. 

Note again, that the- quantities that can be- eemiputeel are- not the- absolute- values of 
Zi(t) but the relative- values i,-(f) — Zj(l). It. follows that withemt. loss of generality, we can 
set (£(),?/(),£o) at the- origin of Mu- coordinate system, that is (:««.;</o,i-») - (0,0,0). At 
every instant, the remaining n — t points are- them given relative to this first point. Thus 
the number of hide-pe-nele-nt variable's in the meiele-1 is n — t. Ini rexhicing the- folleiwing 
neitatiemal simpliiieatiem: 
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a tJ (xi - Xj)(.r, - x.j) I (?/, - VjUih - ?//), 
Eq. (3.5) can then he rewrit ten a.s follows: 


(3.(5) 


n 2 n 1 n l 

/;,(/.) = f ( z * ~ z j )(* ~ *i))* + 5Z("ty + -v^) 2 - (3.7) 

i I J •» j J 

As described above, we are looking for the that, minimize D,.(t,). The necessary 

condition for such a minimization is that the partial derivatives of D,-(t) with respect to 
the i,;(t) are zero: 


01 ),. 


0 . 


(3.8) 


In sect ion 3.5 wo show t hat Eqs. (3.8) represent a set of n 1 linear equations with n — 1 
variables i,(t), which generally have a unique solution. It follows that, because D c (t) > 0, 
thus being bounded from below, the i t (/.) that satisfy Eqs. (3.8) also minimize D e (t), so 
that, these equations generally represent a sufficient condition for the minimization. 

The z.i(t) that satisfy Eqs. (3.8) are expressed in terms of Zi{t.), thus representing 
a system of n - i differential equations with n — 1 variables. This system, however, is 
generally difficult, to solve explicitly, as it is nonlinear. The only solutions that we were 
able to verify by straight substitution are the true motion of the object and its depth 
reversal (which are equivalent, under orthographic projection). In the next section we 
present a stability analysis of the true rigid motion solution; that is, we examine the 
stability of the algorithm when its solution is close to the true solution. (The computer 
simulations address the 1 full convergence behavior of the algorithm.) We concentrate on 
rigid motion, because we feel that nonrigid transformations may disrupt the stability 
of the solutions, introducing instabilities due to peculiar typos of motion. For example, 
the internal model may be too sluggish in adjusting itself when fast changes of structure 
occur. By considering rigid motion in this analysis, we address more fundamental aspects 
of the theory itself. 


3.2.2 Stability Analysis of the Continuous Formulation 

The idea of a stability analysis can be stated as follows. Suppose that at a given instant to, 
the computed 3 D model is very close to the true solution, that is z l (t l j) = ^i(to) + c i{h) 
where z,:(t) is the true depth value at point i. Because the system is perturbed at to, we 
expect it in general to be perturbed for every t > /, (i , that is Zi{t) = Zi(t) + The 

system is said to be asymptotically stable if the following is true: 

lim (i{t) — 0, (3.9) 

t > OO 

for every i. If the (i(t) remain bounded, but do not converge to zero, the system is 
defined to be weakly stable. If, however, lim t >oc <,;(/,) ~ oo for some i, the system is said 
to be unstable. 
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If we suppose Mull < , (/.) is small enough for every I > t. o. t hen we can make I he 
following first order a,pproximat ion of' Eqs. (3.8) in terms of <,(/.) and <',(/) (tin- arguments, 
/,, have been removed for simplicity): 


/ D-D,. D-D,. . \ 


(3-iO) 


A derivation of Eqs. (3.10) is given in Appendix B. Eqs. (3.10) can be readily solved for 
<j(t) in terms of <y(/,). For nofational simplicity, we write this solution in matrix form: 


' i)-D c ' 

i 

i- 

c 



. J 


t — A< 


(3.11) 


w1hto£ ™ ((j,o.,( fl i) ", \() 2 D c / ()z{()zj] is tlio matrix whos<M'l('in('ui at row i and column 
j is d 2 D c / <) i;() Zj and \i) 2 D c /<)z t t)zj\ 1 is the inverse of the mat rix whose element at 
row i and column j is <) 2 ]) c /i)z.,<)zj , provided it exists. The matrix A is defined to be 
D,,/dii<)zj] 1 [D^Dc/OziDzj], and is evaluated at tin' I rue solution itself at every 
instant. 


Using the Eqs. (3.11) to study the stability of the system, it is not possible to prove 
that the system is unstable. To do so requires a proof that, r is unbounded as t —* oo, 
but in this case the first order approximation used to derive Eqs. (3.10) no longer holds. 
It is still possible to determine, however, whether or not the system is asymptotically 
stable. Indeed, we will show in this section, that for many types of motion, the continuous 
formulation of the incremental rigidity scheme is not asymptot ically stable. The results of 
computer simulations provide evidence that the formulation is, in general, not unstable, 
thus being weakly stable. 

The Eqs. (3.10) represent a system of ordinary linear differential equations, which 
are used in conjunction with Eq. (3.7) where D c (t) is defined. Note therefore, that A 
has time dependent components. How do the different types of motion determine A? A 
generalized motion of a rigid body can be described instantaneously as a rotation about 
a fixed axis in space plus a translation of the origin of the coordinates. We noted before 
that under orthographic projection, all that can be determined are the relative depths of 
the points and not their absolute values. Also, in the same case, the only relevant data 
to the problem of recovering structure from motion are the relative image coordinates 
(terms ci,y in Eq. (3.7)). Translation does not change any relative distances, as if. changes 
the positions of all points by the same amount. Thus only rotations need to be considered 
in this problem. 

As mentioned earlier, the matrix A in Eq. (3.11) is in general time dependent. Thus 
this system cannot be solved by the standard method of characteristic values, available to 
systems with constant coefficients. Furthermore, the general rotation can have variable 
angular velocity, and the axis of rotation can change over time, making Eq. (3.1.1) very 
difficult, to integrate analytically. 

It is not necessary, however, to solve Eq. (3.11) in order to reach certain conclusions 
regarding the uouconvergence of <_ to 0. For example, let f ■ (t) be n --1 arbitrary solutions 
of Eq. (3.11) with initial conditions <; ■(()) specified. Then the n 1 x n — i matrix E(/,) = 






22 


[< I ,(/.)] satisfies 1 he Liouviile Jacobi formula. (Yakubovich and SIarzhinskii, 

1975): 

/h/(E(/)) I)cl.(E(0))cxp(^^ Tr(A(t'))<lt'^J . (3.12) 

11, follows Dial a, necessary condition for Jim* >00 E(<) 0, that is, for the system to be 

.asymptotically stable, is: 

/• OO 

/ Tr(A{L'))dt' = -oo. (3.13) 

./() 

Wo found that condition 3.13 does not hold for many interesting typos of motion. 
For example, for any three point configuration, rotating with arbitrary angular velocity 
about, an axis perpendicular to the viewing axis, we can derive the following relationship, 
by using Eqs. (3.7) and 3.1 1 (details of the derivation are given in Appendix C): 


(2i'> - 

Tr(A(t)) - 


z { )z 2 + (2ij - z 2 )zi 
z\ - iji '2 + z\ 


(3.14) 


where zi, Z 2 , zi, %2 denote the z coordinate’s and velocities evaluated at time t. This trace 
is periodic and has zero integral over the cycle, for rotations of fixed angular velocities 
u>: 


p2ir/oj 

/ Tr{A{t))dt = 0, (3.15) 

Jo 

a result that contradicts Eq. (3.12). Thus, for any three point configuration have this 
type of motion, the internal 3 D model computed by the continuous formulation will not 
converge to the true structure. 

The fact that this result is derived for three point structures is nontrivial. In fact, it 
is shown in the analysis of the discrete formulation, later in this section, that converging 
models can be constructed for three point configurations under certain conditions. A 
detailed analysis for structures with a higher number of points is in general very cum¬ 
bersome, but we verified that Eq. (3.15) holds for special four point configurations, 
such as those that when viewed along the rotation axis, appear as squares, rectangles or 
trapezoids. 

We also note that having the rotation axis perpendicular to the viewing axis is the 
best situation, as far as convergence is concerned, because depth motion is lost as the 
rotation axis is slanted towards the viewing axis (that is, away from the imago plane). 
In particular, note that in the limiting case, where the rotation axis is parallel to the 
viewing axis, no motion in depth occurs whatsoever, and so the structure cannot be 
recovered from relative motion. 

Is the instability of the continuous formulation due to properties of the projections of 
the object, at particular positions in its revolution, which somehow convey less information 
to the recovery of st ruct ure front motion (for example, when many of the points belonging 
to an object, have little motion in depth)? In order to check this hypothesis, we explored 
the convergence of the internal model when the object, was oscillated back and forth 
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over a small angular extent,, around different positions. It was found again, by (becking 
the validity of Eq. (3.15), that regardless of position, the internal model is unable t.o 
converge to (lie I,rue structure, t.lms giving a negative answer to the question raised. 'Phis 
suggested to us that the reason underlying (In' instability of the continuous formulation 
may bo based on the smallness of the angular displacement between consecutive frames, 
and not on the structure of’the frames themselves. Pursuing this idea further, the next 
sections address the convergence behavior of t he discrete (and more general) formulation 
of the incremental rigidity scheme, focusing on the effects of the size of the angular 
displacements between frames. 


3.2.3 The Discrete Formulation 


The continuous formulation was derived from the approximation of Eqs. (3.3), as the 
angular displacements between the consecutive frames used in the structure from motion 
process become small, and it was expressed in Eq. (3.5). If we do not make such an 
approximat ion we have a more general formulation, which is however discrete. By being 
more general, this formulation allows a study of the incremental rigidity scheme as a 
function of the angular displacement, between frames, and enables us to verify whether 
the total instability found in the continuous formulation is due to an “artifact” of the 
approximation made, namely an infinitesimal angular displacement bet ween frames, or 
rather is due to a continuous process in which the system becomes less stable with smaller 
angles and eventually becomes unstable when the displacement is small enough. 

In order to stress the discrete nature of the general formulation, we first rewrite Eq. 
(3.3), using as independent variable! t = kr, where k is an integer and r a fixed interval: 

D d (kr) = Y2 [fai(fcr) - zy(fcr )) 2 -f (?y;(fcr) - Vj{kr )) 2 + {zi(kr) - Zj{kr )) 2 
id 

~ ( x { {(k -I- l)r) - Xj-((k + l)r)) 2 - (?y;((fc + l)r) - yj{{k + i)r)) 2 ^' 16 ^ 

- (zi((k + l)r) - Zj((k + l)r)) 2 ] 2 . 

Note that as in the continuous formulation, the quantities of interest for the theory are 
relative rather than absolute, i.e. (zi((k + l)r) — Zj((k + l)r)). Thus again, without loss 
of generality we can sot (:r ( ), ?y ( ), 2 (J ) = (0,0,0), and study the behavior of the other n — 1 
points relative' to it. For notational simplicity, let: 


bij(k T ) =(xi((k -I- l)r) - Xj((k b l)r)) 2 + (y-((k + l)r) - yj({k + l)r)) 2 
- (xifArr) - x.j(kT )) 2 - (yy,:(fcr) - ?/j(A:r)) 2 . 


We noy; rewrite Eq. (3.10), using this notational simplification and without explicitly 
specifying r: 




24 


n 2 n J 

u,(*) EEfM*) 1 (-.■(* f 0 I I)) 3 (-.(*) M*)) 3 ) 5 

* 1 i -i 

n 1 

i £>>,•(*) ' * 3 (*' l ) -*,w. 

j j 


(3-18) 


In I,ho prosont. formulation we arc looking for Mio Zj(k I .1) Ihaf minimize D,i(k). As 
in the continuous ease, the necessary and suflicient condition for this minimization is: 


iU) d 


sr (), 


(3-19) 


<>~.,(k t 1) 

which represents a set of n I nonlinear equations implicitly relating the Z{(k) and 
Zi(k I l). The only instance in which we have been able to derive a full analytic solution 
for these equations is in the two point cast'. This case is a degenerate one, because for 
two points the 3 D structure is not determined uniquely by any number of views. It is 
important to discuss this case, however, as it provides a relevant insight into the theory, 
which is discussed further in section 3.2.5. Eqs. (3.18) and (3.19) in this case reduce to: 


M*) + 4( k + 1) - z\ (h)] z x (k + 1) = 0. (3.20) 

The only nontrivial solution for this equation occurs when the term inside the brackets is 
zero. But. using the definition of I) ( n, this condition essentially implies that the distances 
between the two points is constant. In other words, the distance 1 between the points in 
the initial internal model at k — 0 will tend to remain the same 1 throughout the entire 
motion even if it is wrong. The model will only correct itself if its initial guess for the 
distance between the points is small. In that case, the internal model will be forced 
to change because the distance between the points in the projected plane will become 
larger than the distance in the model. In this case, numerical calculations show that 
the model expands, and even has a small overshoot such that the distance between the 
points becomes a little larger than the true distance, and then stays at this condition 
forever. 


3.2.4 Stability Analysis for the Discrete Formulation 

The stability analysis for the discrete formulation begins with the same general argu¬ 
ment as for the continuous case. Suppose that at k — 0, the model is very close to the 
true solution, i.e. z,-(0) = £ t -(0) I-1,(0), where £,•(&) is the true depth value. We exam¬ 
ine whether or not the perturbation at later times converges to zero, that is, whether 
lim/ c ,oo u{k) = 9. If this is the case, and the ( 7 :(k) always remain small, then we can 
write an equation similar to IOq. (3.11) for the discrete formulation: 


i)*D d 

l 

i^Da 

i)z. L (k -\~ L)dz.j{k -fl) 


<)z.i(k -f- 1 )dz.j{k) 


c(k) = B(k)<(k). 


(3.21) 


1 (k + 1) = - 
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Again, by the same arguments given in the continuous cases, we cannot use Eq. (,‘5.21) to 
prove the system to be unstable, but only verify whether its convergence is asymptotic. 

In this analysis we concentrate on three point configurations, because as in the con¬ 
tinuous formula! .ion, an analytic st inly with a larger number of points is too cumbersome. 
In the results of computer simulations, however, we found that the results derived here 
seem to generalize readily to configurat ions with a higher number of points. We found in 
all the cases studied, that, if the viewed object is rotating with constant angular velocity 
around an axis perpendicular to the viewing axis, then the discrete formulation yields 
a converging internal model, but that the rate of convergence' depends on the angular 
displacements between consecutively viewed frames. 

In order to define the concept of rale of convergence, we first point out that recursive 
use of Eq. (,1.21) shows that < (k) can be readily computed from <_(()): 

k 1 

L (k) =-- (J1 B(j)k(0), (3.22) 

i o 

and thus the convergence of < (h) to 0 depends only on this matrix multiplication, rather 
than on the perturbations themselves. Furthermore, if B(j) is cyclic (as it is for rotations 
in which the ratio between the angular displacement and 2ir is rational), with cycle m, 
convergence depends only on the multiplications of the matrices over the cycle, which 
we denote: 


rn 1 

B m = n B (i)- ( s - 23 ) 

j 0 

Interestingly, the behavior of Ii ] m bears a strong similarity with geometric sequences. 
Let us define the spectral radius of B m . p( B m ), as the maximal modulus of the charac¬ 
teristic values of B m . Then it is possible to show (see Appendix D) that: 

P y (B m ) = p( B{ n ). (3.24) 

From this result, one can understand the following important conclusion (see, for exam¬ 
ple, Varga, 19G2): 

bm B l 0 f>( B m ) < 1. (3.25) 

3 ->oo 

In other words the necessary and sufficient condition for the' internal model to be asymp¬ 
totically stable around the true solution, is that the spectral radius of the rotation matrix, 
B m , is less than 1. 

Furthermore, from Eq. (3.24). one sees that the spectral radius corresponds to a 
“time constant.” by which we can measure the velocity of convergence of the internal 
model. In fact, as p becomes closer to 1, more revolutions are necessary to obtain 
the same amount of convergence of the system, and the spectral radius of B 3 n , declines 
exponentially with j. This exponential decline, however, holds only for the case of small 
perturbations, from which Eq. (3.22) was derived, and describes only a general trend, as 
in line detail, the cyclic matrix B m is composed of a multiplication of parthil matrices 
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(hat may impost* a “noise” in tin* decline. This noise can lx* set'll in (he simulations 
shown in sect ion 3.1. 

The dependence of (he rate of convergence of (lie internal model on the angular 
displacement, between viewed frames, <1, and its consequences, art* better illustrated by 
some examples. Fig. 8 displays p as a function of d for the cast* in which I,he viewed 
object appears as an equilateral triangle, when viewed along the rotation axis. 



--- - Angular Displacement 

Figure 8. Graph of the spectral radius p as a function of d for the case in which the viewed 
object appears as an equilateral triangle, when viewed along the rotation axis. 


The first result of interest, in this graph is that p < 1 for every angular displacement 
d such that 0° < d < 90°. This means that for almost every displacement used by the 
algorithm, the true solution (up to a depth reversal) is an asymptotically stable one. This 
result held for every configuration tried, provided that the rotation axis was perpendicular 
to the viewing axis. Thus, in spite of the fact that the true solution is asymptotically 
stable in the cast* of the discrete formulation, in the limit of arbitrarily closed frames, 
that is, in the continuous formulation, the true solution is not asymptotically stable. As 
can be seen in the figure, the spectral radius tends to .1 as d 0. Thus, in the limit, 
the condition expressed in Eq. (3.25) for the convergence of Bj u to 0, is not fulfilled. 
This confirms onr previous results regarding the lack of stability for the continuous 
formulation. 

In addition, from the limit of 1 at small displacements, p decreases with d, showing a 
steady improvement in the rate of convergence of the internal model to the true solution 
with displacements up to GO". This observation was also made by TTlIman (1984). The 
behavior of p with small displacements is examined more closely in the next section. In 
the rest of this section we discuss the affect, of largo angular displacements on the spectral 
radius. 

For large displacements, the rate of convergence deteriorates as can be seen in Fig. 
8, approaching the value of 1 for d — 99°. The spectral radius is J for displacements 
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of 90' because only two different views of the object are available under orthographic 
projection, which is too little information for recovering structure from motion based on 
the assumption of rigidity. This result appears to be independent of other properties of 
the moving points. 

This behavior of p with angular displacement, can be related to the convergence rate 
observed in computer simulations. This is illustrated in Fig. 9. We began with the 
set. of three points, whose true coordinates (when projected onto the ,Y — Z plane) lie 
around an equilateral triangle, as shown in Fig. 9a. The y coordinate's of the' points were 
nonzero but small. These points were then rotated around the vertical axis, and the 
error measure shown in Eq. (3.1) was computed, as a function of time. Fig. 9b shows 
an error plot for the case in which the angular displacement, between frames was 20 , 
and the points were rotated for a total of 10 revolut ions. The corresponding graph, with 
error plotted on a log scale, is shown in Fig. 9c. It can be seen that on a log scale, the 
long term convergence of the computed 3 D model is approximately linear. A measure 
of convergence rate can be given by the slope of the line that best fits this data, in a least 
squares sense. In Fig. 9d, the slope of tins line is plotted as a function of the angular 
displacement, between frames (the three point configuration was rotated for a total of it) 
revolutions in each case). It can be seen that this graph is qualitatively similar to the 
graph of p shown in Fig. 8. 

3.2.5 Convergence Under Small Angular Displacements and Instability of 
the Continuous Formulation 

The question addressed in (.his section is, why does the rate of convergence fall with 
decreasing angular displacement between consecutive frames, eventually leading to in¬ 
stability for the case of the continuous formulation? One would certainly expect poorer 
corrections for the perturbations if smaller displacements are used, because less resolu¬ 
tion is available in the image data. That is, the signal to noise ratio is substantially 
reduced. Many more views, however, are seen in the case of small angular displace¬ 
ments, which might trade off with the deterioration in the data resolution. One must 
keep in mind that there are two types of noise that can influence the performance of the 
incremental rigidity scheme: (1) the noise from the image data, and (2) the difference 
between the model and the true solution. In the computer simulations presented ear¬ 
lier, the image measurements were computed analytically, and hence had little error, so 
that the primary source of noise was the discrepancy between the computed model and 
true structure. This source of noise contributes to the deterioration of the solution with 
smaller angular displacements between frames. From the results in the previous section, 
we can conclude that in some sense, with decreasing displacement, the deterioration due 
to poorer corrections occurs at a faster rate than the improvement with the number of 
views. We explore this phenomenon in this section. 

Formally, we can express the above question as, why does the spectral radius p 
of the rotation matrix B„, increase with decreasing angular displacement, d, eventually 
tending to 1 when d — > 0? This is certainly a characteristic of tin' curve of p vs. d 
curve if d is small enough (see Fig. 8). The matrix B m results from the multiplication 




onto the X — Z plane, (b) The error in the computed 3 D model, as a finic.tion of time, for 20 
displacements between frames., (c) The error in tlie 3 D model, plotted on a log scale, (d) Tk 
rate of convergence of the solution, as a function of the angular displacement between frames. 
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of (In' matrices corresponding to partial displacements (Eq. (,‘5.2.'!)). One would like to 
define souk 1 quantity related to these partial displacement matrices, which expresses (.lie 
deterioration of their information content, and to study the connection of this quantity 
to the general spectral radius as a function of d. This quantity and its relationship to 
the spectral radius were found by numerical experimentation. The intuition underlying 
these* observat ions is what follows. 

(liven throe point configurations and fixed angular displacements, an inspection 
of the partial displacement, matrices shows that they all have* exactly the same two 
characteristic values, c j and r 2 , which are real, and 0 < eq < 1 and 1 < e 2 . If these 
characteristic values belonged to a diagonal matrix, then the effect of such a matrix would 
be to change the perturbation vector < (<j.< 2 ) 7 into (<q<q,c 2 c 2 ). hi general, however, 

the vector < can point in any arbitrary direction, and the effects of a diagonal matrix 
on the length of < would be different, from direction to direction. Thus a meaningful 
quantity belonging to these matrices must be some average* of tin- effect of the matrix 
over all directions of the perturbation vector. An important example is tin* mean length 
squared of the vector resulting from the application of the matrix to the general unitary 
vectors, (cosq, smq). This quantity is tin* following: 

1 r 2,r c 2 4- r 2 

i- J o (<jro, 2 ( 7 ) + c!«m 3 ( 7 ))rf7 = (3-20) 

The partial displacement matrices are in general not diagonal, but. if d is small 
enough they can be closely approximated as orthogonally similar to the diagonal matrix 
having eq and e 2 as its elements. If we begin with an arbitrary perturbation vector, 
these matrices will rotate the vectors and transform the resultants as described above, 
and then repeat, this transformation until a cycle* is completed. Thus for small enough 
displacements, one can expect the diagonal matrices to work on perturbation vectors of 
almost all directions, thus having an average effect as described in the last paragraph. 
Generally this mean effect always represents a decrease in the perturbation vector, which 
is applied 2k/ d times during a cycle of revolution of the object. Thus for small enough 
displacements one ran expect a connection between eq, rq, p and d given by an equation 
of the following form: 


p{d) = 


c\{d) + 4(d) ^ 


2rr/d 


(3.27) 


Eq. (3.27) closely matches the behavior of the spectral radius, for the equilateral 
triangle case (Fig. 8) for d < 20°. For other configurations, relationships very similar to 
this one wore observed to apply. The importance of Eq. (3.27) is that it relates a quantity 
that can be computed in any partial displacement matrix, namely (cj(d) +c 2 ( ( 7))/2 (the 
deterioration factor), to the spectral radius of the full rotation matrix, and it enables us 
to reduce the problem of why rite continuous formulation is unstable to the question of 
how fast does the deterioration factor tend to 1 as d —* 0. In fact one can see from Eq. 
(3.27) that a necessary condition for liny ,o p{d) — I is that: 



30 


lilt! 
d > 1 ) 


4 (d) 1 <‘?>(<l) 


2 


l. 


(3.28) 


This con (lit ion, however, is not sullicie-nt. For example, if the deterioration factor 
approaches 1 linearly with d, i.<'. (c((r/) I (, |(d))/2 m I ad for small tl, then because of 
the power 2n/d in ICq. (3.27) w<' would have lima .» f>(d) <: < 1. We expect then, 

that the deterioration factor approaches 1 faster than linear. This is indeed the case, as 
can be seen by expanding ( c\(d ) I 4(d))/ 2 into the first few powers of its Taylor series. 
For any three point configuration, this expansion yields: 


i 4(d) 


l - /hr. 


(3.29) 


This quadratic approach of the deteriorat ion to 1 is fast enough to ensure t hat lima >o “ 
1, as can be seen by substituting Ee). (3.29) into Eq. (3.27) and taking the limit. 

Therefore, we confirm that the instability of the continuous formulation is due to 
the fact that the deterioration in the data resolving power as a result of the small angular 
displacements between frames is faster ((juadrat.ic with the inverse of the displacement) 
than the increase of informat ion due to the increase in the number of frames (linear with 
the inverse- of the displacement). This analysis provides theoretical support for the ob¬ 
servation that for the case of the discrete formulation of the incremental rigidity scheme, 
the rate of convergence of the internal model to the t rue solution diminishes as the an¬ 
gular displacements between frames decreases. In the limit of very small displacements, 
that is, in the case of the continuous formulation, this deterioration leads to a lack of 
asymptotic convergence of the model to the true structure. 


4. Perspective Formulations of the Incremental Rigidity Scheme 

In this section, we first present discrete and continuous formulations of the incremental 
rigidity scheme that use perspective projection of the scene onto the image plane. Wc 

then present a few 4 * * 7 * * 10 example's of the results of computer simulations, in which the per¬ 

spective formulations were applied to sequences of points undergoing both pure rotation 
about a single axis in space and pure translation through space. The behavior of the 
perspective- formulation of the incremental rigidity scheme is more complex than that 
of the orthographic formulation. We- found that if the absolute position of an object 
in space is known throughout the motion of the object, then the perspective formula¬ 

tion performs well, similar to the orthographic formulation. If the absolute position of 
the object is unknown, the incremental rigidity scheme in general does not derive the 
correct 3 I) structure. The results of computer simulations revealed a degradation in 
performance with smaller angular and spatial displacements between frames, but this 
de-gradation was somewhat more severe than in the orthographic formulation. In the 
limit of the continuous formulation, the solution is no longer stable. 

In both the discrete and continuous formulations, we assume- that the positive z axis 
points in the- dire-ction of the optical axis, with (lie- image plane at 2 = 1, as shown in Fig. 

10. For notation;),1 conve-nience, we let (£,(/,), ;(/,•(/), 2 ,(/,)) repre-seut the true coordinates 





4.1 The Discrete Formulation 


In the discrete case, we compute the depth values Zi(t') that minimize the measure of 
rigidity, D^(t,t') given by the following (for notational simplicity, (x{, i/i, 2 ,) refer to the 
coordinates of the point s at time t in the computed model, while refer to the 

coordinates of the points at a later time t '): 


Ai(M') = £(M<) - MO ) 2 

».y 

= [(( X * “ X i ) 2 + (Vi ~ yj ) 2 + ( Z i - z j)*) 5 

- (( j n- - x 'j ) 2 + {y'i - y'jf )- (s'i- z 'j )“) v 


(4.2) 
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'['his measure is (he same as (hat used in the orthographic formulation (i.e., it is derived 
from ICq. (2.1) with the l'L(l-) factor in the denomonator omitted, as shown in Eq. (Al) 
of Appendix A). 

As before, we assume that we have a current model (:/:;(/,),?/,;(/,), ;:,:(£)). Rather than 
assuming that. Xi(t') and are known, howi-ver, we assume that, only «,•(/,') and ?;,,;(£') 

are known (they are derived from the true spare coordinates at time The coordinates 
of the points in space at time for the computed model are then given implicitly by: 


Xi(i !') = 

Hi(t') -- ».■(*')*(*')• (4.3) 

By substituting Eqs. (4.3) into Eq. (4.2), D,i(tj') can be <*xpressed in terms of the 
coordinates of the points as follows ((-«',u') demote' the image coordinates at time £'): 


D <i{ t, I') = (( x i - x j ) 2 + (»t - Vj ) 2 + («i - Zj) 2 ) 2 






U 


3 3 


(4.4) 


) 2 + WA - v 'jz'j ) 2 + U - -}) 2 ) 2 


After the values Zi{t') that minimize Dj{t, t 1 ) are computed, Eqs. (4.3) are used to derive 
the space coordinates and yi(t'). The new 3 D coordinates, (z;(f'),•(/;(£'), Zi(t')), 

then become the current 3 D model, a now frame in the sequence is registered, and the 
process repeats itself. 

In the case of perspective projection, 3 D structure can be recovered from relative 
motion only up t o a multiplicative scale factor. This is clear by inspection of Eqs. (4.1); 
a constant scaling of the space coordinates (x t -(<), 2/«(0> *i(0) ( h»'s not change the image 
coordinates (w,:(t),w t (<)). In most of the computer simulations, we assumed that the 
overall scale factor is known; this is equivalent to assuming that the absolute position of 
the object in space is known. In the simulations, the absolute coordinates of one of the 
n points in motion was supplied to the algorithm and the coordinates of the remaining 
n — 1 points were computed relative to the known point (similar to the simulations of 
the orthographic formulation). If the absolute position of the object in space is not 
known, then the scale of the computed 3 D model depends on the choice of the initial 
configuration. For example, suppose that we assume that the set of moving points 
initially lies on a plane that is parallel to the image plane. The scale of the computed 
3 D model will then depend on where this initial plane is placed in depth. That is, if 
the plane is positioned twice as far from the image plane, the computed 3 D model will 
essentially be twice as large. We will briefly mention the results of computer simulations 
in which the absolute position of the object in space is unknown. In these simulations, 
when we began with a flat; initial configuration, we usually placed the initial z coordinates 
of the points at a depth that was approximately the mean of the true z coordinate's of 
the points in space. Further details of the implementation of the discrete' formulation 
can be found in Appewlix A. 
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4.2 The Continuous Formulation 

In the continuous ease, wo compute the z components of velocity dial, minimize the 
measure of rigidity, given by the following (the coordinates and velocities are all 

measured at time t): 


D,(‘) 


£ 


[Xi 


- i j > + ( y« - 7/j )(?>,• ilj) -I {~i -- g j)( j,- - h ; )j' 

(*t - -^) 2 + (?/z - '//.,) 2 4- (2i - Zj) 2 


(4-5) 


This measure is the same as that used in the orthographic cast' (set' Eq. (2.7)). It is 
assumed that a current model of the coordinates of t he points in space* is given, along with 
the known coordinates and velocities of the points in the image plane, «,(/,). ?;,;(/.),«.;(/.), «,(<) 
The velocities of the points in space for the computed model are then given implicitly 
by: 


Xi{t) - Ui(t)Zi(t) + Ui(t)zi(t) 
Vi(t') = Vi{t)zi{l) + Vi{t.)zi(t). 

If we let. Xij — Xi - Xj, y.;.j = ?/; - j/y, Zij — Zi - zy, and i;y 
expressed in terms of t.ht' coordinates of the points as follows: 


(4.0) 

Zi -- Zj, then D c (t) is 


D e (t) = E 




(UiZi 


“J "3 I 


) + Vi:j(viZi + Vi 


VjZj) + Zi. 


1,3 


X?- -+ 

l] 


1/?. + 2?. 
■>IJ l) 


(4.7) 

In other respects, the continuous formulation is similar to the discrete formulation. A 
model of the structure of the moving object is built up by continually taking into account 
new velocity information over an extended time period. After the z t ( t. ) are computed 
by minimizing Eq. (4.7), Eqs. (4.0) are used to derive X{(t) and yi(t), which are then 
ust'd to compute the new 3 D model. Again, because perspective projection is used, the 
structure of the object can only be computed up to a multiplicative scale factor. Further 
details of the implementation of the continuous formulation are presented in Appendix 


A. 


4.3 Computer Simulations 

In this section, we present some results of the computer simulations of the discrete 
and continuous perspective formulations of the incremental rigidity scheme, for a small 
number of examples. We first consider the case in which a set of discrete points is 
translated through space, and then consider the rotation of a set of points about a 
vertical axis. Through simulations of the discrete formulation, we examine the rate of 
convergence of the algorithm and the quality of the solution that, it yields, as a function 
of the size of the angular and spatial displacements between frames. We then observe 
the limiting behavior of the continuous formulation for the case of rotation. Wo also note 
the difference in performance of the perspective formulations when the absolute position 
of the moving points in space is either known or unknown. 
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It should ho stressed that in the case' of orthographic project ion, it is not possible 
to recover 3 1) structure for objects undergoing pure translation, because the relative 
positions of projected points in the image do not change as objects translate. The use of 
perspective projection has the advantage of allowing tho recovery of structure for trans¬ 
lating objects. In the simulations presented here, objects were oscillated back and forth 
in the direction parallel to the image plane. If the objects wore allowed to translate 
in one direction for a large extent, the projected image' of the points would eventually 
become so small that the recovery of structure would become difficult due to a loss of 
numerical accuracy in the measurement, of the changing positions of the points. Oscil¬ 
lating the points allows us to analyze their structure over an extended time period while 
maintaining a relatively large iimage. 

In the' first set of examples, the input consisted of a sot of six points in space. The 
coordinates of flit' points were chosen randomly and the coordinates of one of t he points 
was knowm t.o the* algorithm. The projections of the points onto the X - Y and X — Z 
planes are shown in Figs. Ida and lib, respectively. This set of points was oscillated back 
and forth in a direction parallel to the image' plane. The X - Z projection (birds’ ('ye 
vi('w) of the rightmost and leftmost position of the points in space are shown in Figs, lie 
and lid, respectively. The x coordinate's of the points in Figs, lie and lid are shifted 
by GO units from the initial positions shown in Fig. lib. In the simulations, the initial 
points shown in Fig. lib were first translated to the right in constant discrete steps and 
then translated to the le'ft with the same' discrete steps. A 3 D model was built up over 
several oscillations of the points. The initial model for the structure of the points was 
Hat, with the £ coordinate's placed ne>ar the mean of the true 2 coordinates of the points. 
For the particular set of points used in these examples, Zi(t) = 80, i = 1,..., n, in the 
initial 3 D model. Once the initial Zi(t) are specified, the initial x,:(f) and yi(t) can be 
determined from the image coordinates U{{t) and v t {t). The one point whose' position is 
knowm throughout the oscillation of the points is indicated by an arrow in Fig. 12a. 

Fig. 12 shows a birds’ eye' view of the computed 3 D model after 1 (left column), 
4 (middle column) and 10 (right column) complete oscillations of the points, using three 
different sizes of spatial displacements between frames. The true 3-D structure of the 
points is shown again in Fig. 12a. In the case of Fig. 12b, the size of the spatial 
displacements between frames was GO, so that e'ach complete oscillation consisted of four 
frames (i.e., the points were displaced to the right in one step, returne'd to the central 
position, displaced to the left in one step and then returne'd again to the central position). 
In the case of Fig. 12c, the points were translated in ste'ps of 30. so that a complete 
oscillation of the- points consisted of 8 frames. Finally, in the case of Fig. 12d, the points 
we're translated in steps of 5, so that each full oscillation consisted of 48 frames. In Fig. 
12e, we- show plots of the error in the computed models as a function of time, for 10 
complete oscillations of the points. The e'rror measure used here is the same as that used 
in the orthographic case, shown in Eej. (3.2). The e'rror plots for spatial displacements 
of GO, 30 and 5 are shown superimposed, with the displacements indicated above e'acli 
curve. From the results shown in Fig. 12, it can first he' see'!) that the algorithm denis 
essentially converge to the true structure of the points. Second, the rate of convergence 
eef the computed model to the true solution decrease's with smalleT spatial displacements. 




Figure 11. A wet of six points with random coordinates is shown projected onto (a) the X — Y 
plane and (b) the X — Z plane, (c) The rightmost position of the points during their oscillations 
hack and forth in the direction parallel to the image plane, (d) The leftmost position of the 
points. 

For example, after a single oscillation (leftmost column), the model shown in Fig. 12b 
(displacements of GO) is clearly closer to the true structure of the points shown in Fig. 12a 
than the model shown in Fig. 12d (displacements of 5). Tin' three error plots shown in 
Fig. 12e also have periodic variations superimposed on a steady decline in error, with the 
amplitude of the variations increasing with smaller spatial displacements. Thus, in the 
case of the perspective formulation, it appears that the performance of the incremental 
rigidity scheme degrades with smaller spatial displacements when objects undergo pure 
translation through space. This result is analogous to the observation that in the case of 







Fig. 13 illustrates the behavior of the perspective formulation, for the case of rotation 
of the points about a cent ral vertical axis. The set of six points shown in Fig. 1.1 was 
rotated about a seventh point that was placed at the position: 

(xi{t),yi(t),Zi(t)) = (0,0,80). 

The axis of rotation passed through dus point and was perpendicular to the A' - Z plane. 
The position of this central point was known to the algorithm throughout the motion of 
the points and the positions of the remaining six points were computed relative to this 
central point. As before, the initial configuration was flat, with Z{(t) — 80, i = 1 , ...,n. 

The set of 7 points was rotated in angular steps of 40° and 5°. Figs. 13b and 13c 
show the computed 3 D model after total rotations of 40°, 80°, 120° and 3G0° (indicated 
on the left), for angular displacements between frames of 40" and 5°, respectively. The 
true positions of the points are shown in Fig. 13a. In Fig. 13d, the graphs of the error 
in the computed models are shown as a function of total angular rotation, for three 
full revolutions of the points. It can be seen that the rate of convergence and quality 
of the 3 D structure decreases with the smaller angular displacement between frames. 
Although similar to the case of orthographic projection, the degradation in performance 
with smaller angles was somewhat, more severe in the case of perspective projection. 

Fig. 14 illustrates the behavior of the continuous formulation of the incremental 
rigidity scheme, for the case of perspective projection. The same set of points shown 
in Fig. 11 was rotated about the vertical axis and the 3 D model was computed at 
infinitesimally closely spaced frames, using the instantaneous velocities projected onto 
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tlio imago. As in the previous mint ion example's, a sevent h point was placed at. the center 
of rotation of the point s, and the position of this central point was known to the' algoril Inn 
throughout the rotation of the points. Tu Fig. Ha, we compare the true positions of 
the points (left.) with the best solution (right) obtained over four full revolutions of the 
points. Although the model is quite close to the true structure' at this pe>sitiem e>f the 
points, the solution oscillates significantly over an extended time’ pe'riod, similar to the 
orthographic e-ase. A graph e>f the’ error in the' computed mealed over the' four revolutions 
is shown in Fig. 141 >. The’ arrow marks the* point at which the solution shewn in Fig. Ha 
was obtained. Timm is an initial fast convergence toward the true’ structure of the* points, 
but the’ algorithm them oscillates with high amplitude away fmtn the true* solution. We 
ce>neduele' that, similar to the* orthographic case', the' elire’ct, use’ e>f velocity informalie>n for 
the reveweTy of structure by the’ incremental rigielity sche'ine eloe’s ne>t lead t,o a robust 
solution. 

We> should finally note’ that we 1 also e'X])le)re'el se)me' e'xamjde's in whiedi no information 
about the' absolute position of the* points in space is knewn. In this case 1 , we’ be'gan until 
an initial configuration that, was Hat, anel otherwise provided no furthe’r constraint on 
the positions of the* points. The points we're either oscillated back and forth in the 
direction paraded t,e> the* image' plane, or rotated about a vertical axis. We' found that 
the incremental rigielity scheme weiulel semietiine's eleudve the correct solution tmeler tlie»se 
conditions, but in general, t he' algorithm eliel not ele'rive’ the e-orre'et structure, although the 
computeel solution was e'sse'iitially rigiel. The’ computeel solution was not just, a scaled 
or rotate-el version of the true structure, but actually a different one altoge-the’r. This 
suggests that some additional constraint may be required to derive the correct solution. 

Tei summarize, the’ peTspe-ctive formulation of the incremental rigielity scheme ap¬ 
peal’s to perform wedl when the absolute position of the obje-et. in space is known. This 
formnlatiem them allows the recovery of structure for objects unde-rgoing pure translation 
through space, as wedl as rotation. The results of computer simulatiems reveah’d a degra¬ 
dation in performance with smalle-r angular anel spatial displacements between frames. 
This degradatiem appe'areel to be more se-vere' than in the; case of orthographic projection, 
when the points were; re>tate’el about a vertical axis. When the points we're translate;el, 
the;re was a more gradual elecrease; in the rate of convergemcc and quality of the computed 
3 D model with smaller spatial displacememts. In the limit of the continuous formulation, 
the solutiem was no longer stable. 


5. Summary and Conclusions 

In this paper, we st.uelied and generalized the incremental rigidity scheme for the recovery 
of structure from motion. This algorithm, as first proposed by Oilman (1984), assumed 
the visual input to consist e>f a sequence of frames, each containing a finite number of 
points. The scheme maintained an internal model of the structure of the vi(-wed object, 
which was updated from frame to frame, so as to be consistent with the changing image 
and to be as rigid as possible. This internal model was shewn to correctly converge to 
the true structure of the object, for both rigid and nonrigid motions. 



Figure 13. (a) The true positions of the points after rotations by 40°, 80°, 120° and 3G0°. 
(b) The computed 3 D model after rotation of the points by 40", 80", 120° and 300°, for an 
angular displacement, between frames of 40°. (c) The same computed 3 D models, for an angular 
displacement of size 5°. ((d) on next page.) 




Ill tlu' present paper, we focused on an observation made by Tillman, that the rate 
of convergence of the internal model to the true structure increased with the amount 
of change between consecutive frames. A major part of our analysis focused on objects 
that rotated with constant angular velocity around an axis perpendicular to the viewing 
axis, and whose image was formed through an orthographic projection, at equally spaced 
angular displacements. Tlu- complete dependence of the convergence rate on the size of 
the angular change was obtained for a large variety of objects. 

Ullrnaa’s observation was confirmed for small displacements. Indeed, the conver¬ 
gence rate first increased towards a maximum and then decreased, as a function of the 
angular displacements between the frames. The deterioration of the algorithm for high 
displacements can be understood if we note that the number of frames available for anal¬ 
ysis in a given total number of revolutions of the object, or in other words, the amount 
of information provided to the scheme, is reduct'd in inverse proportion to the angular 
changes between the frames. 

The deterioration of the algorithm for small displacements, however, has a more 
complex and surprizing basis. In this case, although the number of frames used in a 
given amount of rotation is high, the spatial resolution between the points is low (see 
discussion on the beggining of section 3.2.5). Our analysis indicates that the latter effect 
dominates for very small angular jumps between frames, reducing the convergence rate. 
Furthermore, we found that this reduction is such that in the limit of infinitesimally 
small displacements, where only velocity information is used, the internal model does 
not oven have a full convergence to the true structure of the viewed object, but only a 
rough approach to the correct, solution. 

In analogy to linear algebra, one can say that the transformations from frame to 
frame of the internal model in the incremental rigidity scheme become ill conditioned as 





Number of Revolutions 

Figure 14. (a) The beat solution (right) obtained over four revolutions of the points is compared 
with the true positions (left), for the case of the continuous formulation, (b) The error in the 
internal 3 D distances between points in the model, as a function of time. The arrow in (b) 
indicates the instant at which the solution in (a) was obtained. 

the angular displacements decrease, i.o. they become sensitive to noise. This property 
is not unique to the present transformations, but occur in other important situations 
such as numerical differentiation. The noise sensitivity of a differentiation of order rn 
increases as n m when the number of points n in a given interval is increased (Strang, 
1970). 

As mentioned before, the orthographic projection is in genera] not physically valid. 
It was used in the present, paper, because it, allows a deeper analysis of tin' phenomena 
studied. This analysis is further validated because these results are not limited to the 
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cast' of orthographic projection, or to the ease in which the object rotates. Under certain 
condit ions, the incremental rigidity scheme can be generalized successfully to formula¬ 
tions that, use perspective projection, and that in this case, the structure of the viewed 
object is recoverable both for rotation and translation. Also, similar t.o the orthographic 
case, the convergence rate of t he internal model t o the t rue st met tire decreased, if smaller 
spatial changes between consecutive frames were used. 

Inspection of previous results by other researchers suggests that the use of informa¬ 
tion that is local in space and time is insufficient, to allow a robust, recovery of st ructure 
from mot ion. We distinguish between two types of local information that are relevant, to 
the problem: (I) spatial locality, which refers to the use of a small number of points or 
feat ures of the viewed object,, and (2) temporal locality, which is the use of a small num¬ 
ber of views of the object. Early algorithms for recovering structure from motion based 
this recovery on a limited number of points and views, and were able to weaver the struc¬ 
ture' of the viewed object when it was rigid. When an object deformed, some algorithms 
could determine its nonrigidity through inconsistencies of the underlying equations. It. 
followed that these schemes are not robust against noise, as noise has similar effects as 
nonrigidity. 

It. can be shown that motion information that is extended in space, but temporally 
local, can be used to provide stable recovery of structure from motion (Bruss and Horn, 
1983; Negahdaripour and Horn, 1985). Perceptual evidence, however, suggests that 
spatial extension is not necessary for the solution of the structure from motion problem 
by the human visual system (Borjesson and von Ilofsten, 1973; Johansson, 1975; Petersik, 
1989). 

On the contrary, Tillman’s incremental rigidity scheme, which was applied locally in 
space, both in its original formulation and in the present paper, uses an internal model 
that is updated constantly, in order to extend motion information over time, and converge 
to the correct solution. Some perceptual studies suggest that this temporal extension 
may be used by the human visual system (Wallach and O’Connell, 1953; White and 
Mucser, 19G0; Green, 1961). 

We showed in the present paper, however, that an algorithm that overcomes tempo¬ 
ral locality docs not, necessarily provide a robust solution to the structure from-motion 
problem. It seems that in this case, a necessary condition for a stable solution is the 
use of views of the object that are significantly disparate. In the limit, the use of pure 
velocity information allows only a rough solution to the structure from motion problem, 
which is less stable over an extended time period. A somewhat similar conclusion was 
derived by Tillman (.1983), who argued that in a temporally local scheme, the recovery 
of structure from the instantaneous velocity field is impossible under orthographic pro¬ 
jection, and that for perspective projection, the recovery is unstable. We suggest that 
the problems with using local velocity information are too large to be overcome by an 
extension of this information over time. In other words, the inaccuracies introduced by 
a velocity bast'd scheme cannot be corrected by the use of multiple views of the object. 

Our results therefore hint that even in a temporally extended algorithm, a memory 
of sufficiently distinct past views or internal models of the viewed object is necessary for 
a robust recovery of its structure. We point out that the incremental rigidity scheme, in 
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its original formulation and in the present paper, used memory of only one past, internal 
model. All hough in principle, a memory of many past models could In' used by the 
algorithm, it is remarkable that such a minimal amount of past informat ion is sufficient, 
to provide a robust recovery of structure. Still, it may be the cast' that the use of multiple 
past views or internal models for establishing a new model of 3 D structure could increase 
the rate of convergence of the internal model. 

In further support of the above conclusions, we note that in order to recover 3 D 
structure from measured image velocities, it is essential that these velocities at least ap¬ 
proximate the true projected velocities of moving elements in space. In general, however, 
it is difficult to derive real projected velocities from the optical flow pattern on the eye 
or camera. Factors such as changing illumination, speeularities, shadows, and rotation 
of an object surface' relative to a light source', can create patterns e»f e»pt.ie*<il flew that eh) 
not correspond to the 1 mil movement e>f feature's on a physical surface*. The elilficulty e>f 
cennputing eorrect projected ve'lexdlies provieles an aelelitional me»tivat.ie>n for ne>t basing 
the re'cem'ry e>f structure 1 from motiem elire'cfly on ve*locity infe>rmatie>n. 

Our analysis raise's a iminbcr e>f issues re'gareling the recovery of structure from 
motion in the human visual system. First, it is not elear whe'tlee'r the- visual system 
achie'ves a stable solution to the structure from me)tie)u problem, or a rough solution 
such as that provided by a velocity based scheme. The more robust solution may not. be 
esse'iitial if we* eonsieler that e>the'r perceptual cue's, such as binocular stereo, may help 
to improve the equality of the 3 D solution obtainable by the* structure frenn motion 
process. Further psycophysical experiine'nts are reepiired to examine* this question. 

A second pe>int of interest is t.iiat our analysis suggests that if the visual system 
incorporates a robust solution to the structure from motion problem!, it must be able to 
match corresponding points in ve*ry disparate views e>f mewing objects. The displaceunemts 
between corresponeliug points may be* larger than the spatial limits proposed for the 
short-range perceptual process found fen- 2 D motiem (Braeldick, 1973, 1974, 1980; Pantle 
and Picciano, 1976; Pctersik, Hicks anel Pantle, 1978; Petersik and Pantle, 1979). A 
similar conclusion, based on eliffeu-ent grounels, was formulated by Pctersik (1980). He 
exple)reel the* sensations e*lie-ite'el by stroboscopic simulations of a rotating transparent 
sphere filled with ranelennly positioned elots. By manipulating both spatial and temporal 
variables in the simulation, he concluded that corresponding elements in consecutive 
frames can be matched ewer spatial anel temporal distance's that e'xceoel the empirically 
determined limits of the* 2 D slmrt range process. Even in this experiment, however, 
points that reach the periphery e)f the sphere have small elisplaeemcnts fre>m frame to 
frame that may fall spatially inside the short range process. It remains, therefore, to be 
established that these points do not provide all the information used by the subjects to 
sense the continuous rotation and internal volume of the sphere. 

If the human visual system is able to match very distant corresponding points from 
disparate views of a moving object, the question of how we solve this correspondence 
problem becomes important (for a review of the motion correspondence problem, see 
Tillman (1981)). This correspondence could be established (*ither through the repeated 
use of a short range matching process, or through the use of an explicit long range 
matching process. In the first case, the short, range process could provide an essentially 
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continuous tracking of (cat urns in the changing image, and the positions of feat ures could 
them be sampled by a longer range tracking process that feeds disparate positions info 
the structure from motion process (a similar idea has been suggested by Ullman (1981)). 
In the latter case, a correspondence would bo established directly from disparate views 
of moving objects, without- the aid of the intermediate short, range process. Pet.ersik’s 
(1980) results may suggest, that if is possible to solve this long range' correspondence 
problem for re-covering structure. This emre-spemelence preiblem cemlel be- solvoel in con¬ 
junction with the structure frenn meetion process; that is, a ce>rre'sponde'nce' coulel her 
chosen that subsequently gives rise- to the- most rigiel 3 D inte-rpre-tatiem (such a match¬ 
ing process re-epiire-s a gleebal ele-eision proceelure', and may lee- we-11 suite-el to solution by 
paralle-1 se-he-me-s such as Ilopfie-ld’s “ne-uronal ne-ts” (1984)). (liven the- elifficult.y e>f se>lv- 
ing this le>ng range corre-sponde'nce problem in general, lieewa-ve-r, it se-e-ms unlike-ly that 
she»rt range measure-incuts wemlel not lx* use-el vvhe-n tlie-y are- available*. 
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Appendices 


Appendix A 

In this appendix, we present some of the details of the implementation of the incremental 
rigidity scheme used to derive the simulation results presented in sections 3 and 4. 


Orthographic Projection: The Discrete Formulation 


Ullman formulated the incremental rigidity scheme as the computation of depth values 
Zi{V) that minimize the measure of rigidity given by D(t,t'), as shown in Eq. (2.2). 
The analyses presented in this paper mainly consider rigid objects in motion, which are 
compact, in the sense that the internal distances between pairs of points do not differ 
much from one another. In this case, the additional distance factor in the denomonator 
of the measure used by Ullman has little influence on the resulting solution obtained 
by the algorithm. We therefore removed this factor in most, of our simulations, and 
minimized instead the following expression (the distances Uj{t) and are expressed 

in terms of the coordinates of the points, and {x^y^Zi) refer to the model coordinates 
at time t, while (x' z , ?/', 2 ') refer to the model coordinates at time t'): 


f) = [((** ~ x j) 2 + {Vi - Vj) 2 + (h: ~ Zj ) 2 ) 2 

hi 

The Zi(t') that minimize D,i(t,t') satisfy a system of nonlinear equations given by: 


dD t 1 

Ozi(t') 


= 0 , 


i — 1,..., n — 1. 


(A2) 


Rather than solving this system of equations explicitly, we solved them implicitly through 
the use of a steepest descent minimization algorithm for Eq, (Al), based on the gra¬ 
dient, of D,i(t,t') (Luenberger, 1973). The components of the gradient are given by the 
following: 


<)D,i 
()zi{t') 



(A3) 


In the case of orthographic projection, only relative depths can be recovered. In the 
implementation, we placed the initial point (a;u(f),?/u(f),^o(^)) at the origin of the co¬ 
ordinate system; that is, (^ () (t), j/o(t), Zu(t)) — (9,9,0). At every instant, the remaining 
n — 1 points were then given relative to this initial point. Unless otherwise stated, the 
initial configuration of points was flat, that is, ,?,(()) = 0, for 1.... ,n — I. From this 
configuration, the algorithm can move toward two equally likely configurations, one be¬ 
ing the mirror image reflection of the other about, the image' plane. In other words, this 
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configuration is located at a saddlepoint in tin' solution space for the problem (the solu¬ 
tion space is an n dimensional space in which a value for D,i(t.,t') is assigned to every 
possible combination of depth values z(/')), at which the gradient of D l{ [l,L') is zero. To 
make this gradient, nonzero, the initial solution was perturbed slightly, thus causing the 
algorithm to move in a particular direction. 

The steepest descent minimization method from nonlinear programming (see, for ex¬ 
ample, Luenberger, 1973) was used to compute the Zi(t') for each new frame. The current 
depth values Z{(t) were used as the initial solution for Zi(t') to begin the minimization 
for the next, frame. The objective function is given by and its gradient by Eqs. 

(A3). Golden section search was used to perform the one dimensional minimization 
within the steepest descent, method. 


Orthographic. Projection: The Continuous Formulation 

In section 2, we presented a continuous formulation, which requires the computation of 
z components of velocity, Z{(t) that minimize the measure of rigidity given by D c (t), as 
expressed in Eq. (2.3). In the computer simulations, wo used the particular measure of 
overall deviation from rigidity given by: 


De{t) — ) „ 


*>J 


[{Xj - - Xj) + (yj - yj){in - y 3 -) + (zj - Zj){zj - Zj)Y 

(Xi - Xj ) 2 + {yi - ?/y ) 2 + (Zi - Zj ) 2 


(A4) 


The Zi{t) that minimize Dd(t) satisfy a system of linear equations given by: 


dDc 

dzi(t) 


-0, 


These equations are of the form: 


i — 1,..., n — 1. 


(A5) 


dPc _ n V ( X i ~ X j)( X i ~ *j) + (Vi ~ yj)(Vi - Vj) + {Zj ~ Zj)(Zi ~ Zj) , _ v _ q 

d*i V ~ (*<-*y) 2 + (»r-»y) a + U--^) a 1 ‘ ” j) 

(AG) 

In the case of orthographic projection, only relative z components of velocity can be re¬ 
covered. As in the discrete case, we again placed the initial point, (x n (t), yo(t)>Zo(t)) at the 
origin of the coordinate system throughout the entire motion, so that (:co (£), t/ () (i), Zo(t)) — 
(xo(t), yo(t), Z(){t)) = (0,0,0). At every instant the remaining n — I z components of ve¬ 
locity were then given relative to this initial point. Unless otherwise stated, we again 
began with a flat initial configuration in which xq(0) = 0, for I,...,n — 1, which was 
perturbed slightly so that the gradient of D c (t) is initially nonzero. 

At each moment, z l (t) were obtained by solving a system of linear equations, for 
which we used the simple Gauss Seidel relaxation method. The initial condition for the 
relaxation was usually the set of velocities computed in the previous iteration and i, = 0 
for the first iteration. To integrate motion information over an extended time period, we 
then.made use of the following approximations using i t -(l), //,:(/), and Zi{t): 
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*,•(/) f = Xi{t f St) 

?/,(/) I !li{t)St - I «) (/tG) 

-i(0 -f = z i{t f m)- 

Once the Z{[t) w<'r<’ computed, hy minimizing D r {t). tlie’ii the' c ,(t -f <S/) could he- ele’rive’el 
from the- current, mode'l Z{(t). The> new mode'] was them taken as A 4{t f St) = ( Xi(t + 
fit),yi(l~jr bt),Zi(t + fit)), and t]ie< process was continued. The time inte>rval, St, typically 
corresponded to an angular displacemiemt between frames on the orde>r of 0.1 elegre’e's. We' 
also e'xpe’rimemte’el wit 1 1 angular displacements up to two orde>rs of magnitude smalle'r, 
and found the ejiialitative behavior of the' algorithm to be' the same. 


Perspective Projection: The Discrete Formulatiem 

In the case’ of pe’rspe'ctive pre»je’ction, the x and y coorelinates at time V can be expressed 
in terms of the known image coorelinates at time t' an el the ele'pth value's te> be computed, 
Zi(t'), as follows: 


Xi(t') -- Ui(t')zi(t') 

yi(t') = Vi(t')zi(t'). ( A7 ) 

The measure of rigidity to be minimized is given by substituting Eqs. (A7) into Eq. 
(Al): 


Ai(<» 0 = i( Xi ~ + “ ? A) 2 + ( z i ~ z j) 2 ) 




(A8) 


(( u 'i z 'i ~ u 'j z j) 2 + ( v i z i - v 'j z 'j) 2 + (4 - z j 


l\2\ 1 


) 2 ) 


The Zi{t') that minimize Dd{t,t') satisfy a system of nonlinear equations given by: 


dD,i 

dzi{?) 


= 0, 


i — 1 ,..., n — 1 . 


(A9) 


In the computer simulations, these equations were solved implicitly through the use of a 
steepest elescent minimization algorithm for Eep (A8), based on the gradient of Dd{t,t'). 
The components' of the gradient are given by the following: 


dD d 

dzi(t') 




ij 


l'.. 


1 


KK-4 - u 'j z 'j) + 4(44 ~ v'j z 'j) + (4 - z j)] ■ 


(ALO) 


The computed Z{(t ! ) are used to derive the new x.i(t') and yi(t') using Eqs. (A7), which 
then provide the new computed 3 D model. This new model also serves as an initial 
solution to begin the minimization for the next frame. 
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Perspective Projection: The Continuous Formulation 

In the continuous case, the x and y components of velocity of the points in space can he 
express('d in terms of the known image coordinates and velocities as follows: 


x,i{t) (7) ! Ui(t)zi(t) 

!k(t) - v { (l)zi[t) I v,{t)zi(t). 

If we let Xij = Xi - Xj, y i3 = ?/».- y 3 , z i3 - Zi ~ Zj , and % 
expressed in terms of the coordinate's of the points as follows: 


(Alt) 


Zi - Zj. then D c (t) is 


DM) = £ 


• ] 2 


|Xjj ( tl t Z/ I 'll 7 Z{ 'llj UjZj) t Vi j{ViZi I V t VjZj VjZj) I Z{jZjj ] 




jf? | «?. | 


The Zi{t) that minimize D,.(t) satisfy a system of linear equations given by: 

DD r 


= 0 i — I,..., n -1. 

Ozi 


These equations are of the form: 
<3D, 


dzi 


2 ^ b t j {x t xj ) I Vi(yi - yj ) ~! ( z t zj )] (), 


(A 12) 


(A13) 


(A14) 


where the arc given by: 


Xij(UiZi + UjZj - UjZj - UjZj ) f yijjvjZi + VjZj 


J 


) h Z{jZ{ 




4 + 4 + 4 


(A15) 


At each moment, the i,(/,) were obtained by solving the system of linear equations given 
by Eqs. (A14), for which we used the simple Gauss Seidel relaxation method. After the 
Zi(t) were known, the Eqs. (All) were used to derive the x and y components of velocity 
in space, i,(f) and yi(t). To integrate motion information over an extended time period, 
we then made use of the approximations given in Eqs. (AG) to compute a new model 
(xi(t + St),yi(t h St),Zi(t + St)). 


Appendix B 

In this appendix, we explain how Eq. (3.10) was derived. 

Eq 3.8 depends on the image coordinates (data) and on depth coordinates of the 
points Zi(t) and their velocities Zi{t) (computed values). These depth variables are in 
general displaced from the true motion values. This departure can be incorporated in 
Eq. (3.8) as follows: 


(j Is • 

7 x . (^5 ~1 "1“ f 1 5 * • • j Zn — 1 "1” €n 1 ? ? • * • j %n~ 1 -1) 

(J Zi 


(Bl) 
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If the c and < an' small enough, as it is assumed in the stability analysis of section 
3.2.2, then the first elements of tin* Taylor expansion of Eq. (Bl) yields the following 
approximation: 


, OD. . 

0 ()z- ^ 5 ***’ 1 1 %n — 1 ) 

V ' 1 VD u - - - -X 

3i 1 3 


{D2) 


n 1 


NT - ' ,, , . x ; . . 

+ / > fj:. B' Zn \) ( j- 

3 1 3 

But the true rigid solution certainly minimizes Eq. (3.7) as it sets it to 0, thus being a 
solution of Eq. (3.8). It follows that the first, term on the righthand side of Eq. (B2) is 
0, from which Eq. (3.10) is concluded. 


Appendix C 

In this appendix we indicate the method by which equations of the type 3.14 and 3.15 
were derived. 

Equation 3.14 is derived by taking the necessary partial derivatives of D (Eq. (3.7)) 
by Zi and Zi and combining them as indicated in Eq. (3.11). These derivatives are 
evaluated at the true solution itself, that is at Zi(t) and £i(t), and because we only 
considered rigid motions in this paper, we could use the following relationship to simplify 
the results: 


^ij b ( z i z j)( z i z j) —‘ fb (Cl) 

Eq. (Cl) is derived by setting the distances between points i and j as constant and 
taking the temporal derivative of this distance. 

The conclusion of Eq. (3.15) is derived by noting that, in a rigid rotation around 
an axis perpendicular to the viewing axis, with constant angular velocity w, Z{ can be 
written as: 


Zi(t) = CiCos{tjjt + (f>i). ((72) 

Then, by direct substitution of Eq. (C2) in Eq. (3.14) and its integration as indicated 
in Eq. (3.15), we derive the stated result. 

The calculations explained in this appendix, though straightforward, are cumber¬ 
some, and were done by using Maesyrna, a computer system for performing algebraic 
manipulation. 


Appendix D 

In this appendix we show that if C is an n x n matrix over an algebraically closed field 
then: 
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^'(C) = p{&), (Dl) 

where p is the spectral radius. 

Let P be such that: 

P 1 CP = J, (D2) 

where J is in Jordan canonical form. It follows that: 

& = PJ^P 1 . (D3) 

Similar matrices have similar characteristic values and thus similar spectral radii, thus: 

p( C) = p( J), (D4) 

and 

p(&) = p(V). (D5) 

Direct inspection of the matrices J 3 show that they have the same type of triangulariza- 
tion as J (that is, upper or lower triangularization), with the elements in the diagonal 
being the j — th power of those of J. Thus it follows that: 

P(V) = P j { J), (D6) 

which together with Eqs. (C4) and (C5) imply the result stated in Eq. (Dl). 



